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SECTION 1 


SUMMARY 


This is the final report for the NASA Lewis funded grant “Characterization of Reac- 
tion Kinetics in a Porous Electrode” (originally titled “Cyclic Voltammetry in a Porous 
Electrode”). The name change was approved in the third-project year to reflect the more 
general electroanalytical nature of the mathematical methodology which was being devel 
oped. 

Our objective has been to develop a procedure by which electrode reaction kinetic 
constants can be extracted from impedance or linear sweep voltammetry measurements on 
a porous electrode. The mathematical theory to interpret such measurements on planar 
electrodes is well known; however, due to the three-dimensional, distributed reaction dis- 
tribution in a porous electrode, it cannot be rigorously applied to this electrode. By using 
a tube-analog model for the pore geometry, we have developed such a methodology and it 
is reported upon in the following sections. 

Our approach has been to develop these tools using simplified redox kinetics (Butler- 
Volmer) as the example reaction sequence. It should be emphasized, however, that more 
complex reaction expressions can be used within the context of the models with appropriate 
modification of the computer programs. We have attempted to write the FORTRAN 
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substitution of alternative kinetic equations possible. 


Two chemical engineering graduate students have worked on the project. Mr. Andrew 
Viner received his MS degree in 1988 for his thesis “A Model of Steady-State and Impedance 
Measurements in Flooded Porous Electrodes,” and currently Mr. John Weidner is working 
on his PhD dissertation “Linear-Sweep Voltammetry in Porous Nickel Oxide Electrodes” 
with an anticipated graduation date of no later than May, 1991. Because Mr. Weidner’s 
work is still in progress, this “final” report is not complete. In the last section of the 
document, the ongoing work is outlined. Upon completion of Mr. Weidner s dissertation, 
a copy will be forwarded to the grant monitor, Mr. N. Hagedorn. 

The following sections of this report are self-contained documents and are presented 
in the chronological order in which the work described in each was done. Each section 
has appeared as a paper in a peer-reviewed Journal (Sections 2-4) or is currently being 
prepared for submission (Section 5). The financial support of the NASA grant has been 
acknowledged in each publication. A short summary of each follows. 

The ohmic resistance in a pore will distort the shape of a voltammogram. As a first 
attempt to quantitatively understand the phenomena, Section 2 reports upon the ohmic- 
distortion effect for thin-layer cells when a reversible redox reaction occurs. A simple-to- 
apply correlation to account for the shift in the peak potential as a function dimensionless 
solution resistance is presented. The theoretical calculations based on porous electrode 
theory provide the basis for the correlation. 
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The MS thesis of Mr. Viner provided an extensive theoretical foundation to calculate 
the steady-state and transient response of a porous electrode. Two-dimensional, analytical 
solutions to the diffusion and Laplace equation in the tube geometry were developed, and 
we demonstrated how the two may be used in a numerical calculation to determine current- 
voltage behavior when coupled through the reaction-kinetic expression. Section 3 reports 
upon a specific electroanalytic measurement, the open-circuit impedance, and illustrates 
how the interaction of ohmic, kinetic, and mass-transfer resistances influence the results. 
With this methodology available to separate the various resistances, kinetic constants can 
be extracted from the measured impedance. 

During the course of developing the mathematics to predict the voltammogram for a 
porous electrode, we discovered a minor limitation of the classic Nicholson- Shain treatment 
for voltammetry at a planar electrode. In their work, the assumption was made that the 
reactant only is present at the onset of the sweep. From our mathematical development, 
we determined that if both reactant and product are initially present at comparable con- 
centrations, the peak potential and current for a reversible redox reaction are significantly 
changed from that predicted by Nicholson and Shain. In Section 4, we report upon this 
work which is of importance to the electroanalytical community. The magnitude of the 
effect is a function of the ratio of the initial reactant-to-product concentration ratio times 
the square root of the ratio of the two diffusivities, and a correlation for the peak potential 
and current is presented as a function of this variable. 
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The peak potential and current for linear sweep voltammetry in a pore when ohmic, 
kinetic, and mass-transfer resistances may all be of importance is presented in Section 5. 
Again, correlations are presented, based upon the theoretical calculations, which enable 
the determination of the kinetic constants. This Section is nearly ready for submission for 
peer review except a few additional calculations must be made, as appropriately noted in 
the Section. 
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ABSTRACT 


A continuum- model approach, analogous to porous e.ectrode theory, has been 
applied to a thin-layer cel. ol rectangular and cylindrical geometry. A reversible redox couple 
is assumed, and the .oca, reaction current density is related to the potential through the 
(ormuia of Huhbard and Anson for a uniform!., accessible thin-layer ceil. The placement of 
the reference e.ectrode is aiso accounted for in the analysis. Primary emphasis is placed on 
the effect of the solution-phase ohmic potential drop on the voltammogram characterist.cs. 
Correlation Rations for the peak-potentia, dispiacement from * and the peak current are 
presented in terms of two dimensionless parameters. 
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Unman, Pons, and Cassidy |1] discussed the effects of solution-phase resistance on 
lhe sh a P e of the voltammogram o, a reversible redox couple in a thin-layer cell. In this cell, the 

and a reference electrode was placed between the leading edge of the worhing electrode and the 
counterelectrode. Hubbard and Anson’s [2] formula for the current in a uniformly acceasible, 
thin-layer cell was used to represent the surface impedance in an equivalent network model. 
Both a rectangular and a circular geometry for the working electrode were considered. In 
both cases, the peak potential of the reversible vo.tammogram was shown to be displaced 
from E°' — where it otherwise would be located-by solution-resistance effects. The circular 


5-electrode model [3] of the cell 


was less distorted. 

In this communication, a one dimensional, porous- 
will be used as an alternative to the network model of Hinman el al. The current voltage 
behavior was found, as expected, to be identical in the two approaches. However, in formu- 
, ati „g the differential equation of the present work, two dimensionless parameters are shown 
to govern the behavior of the voltammogram. We present correlations which show how the 
peak-potential displacement from E°‘ and peak current are dependent upon these two param- 
eters. With such relationships available, a thin-layer cell may be, for example, des.gned w.th 

a known deviation from ideal behavior. 
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MODEL FORMULATION 


Consider the rectangular thin-layer cell shown in Figure la with the reference 
electrode located at x = -x RE . Because the cell gap is much smaller than the electrode length, 
the potential distribution in the electrolyte is assumed to be one dimensional and dependent 
only upon x. A differential conservation of charge equation may be written as 



( 1 ) 


where h is the gap height, i 5 is the current density in the solution phase, and i n is the reaction 
current density— with a cathodic reaction current being positive. Subject to the assumptions 
and limitations discussed by Hinman ct al, the local reaction current density is calculable 
from Hubbard and Anson’s formula for the current as 

( nFfuC° ox h e*p\m{E ~ E°'j\ (2) 

In * “ ' {l + exp^E -^')]} 2 

Because the potential is position dependent, so must be i n according to equation 2. If ohm’s 
law is used to relate i, to the solution potential gradient (i, = kV E where E is the potential 
difference between the working and a reference electrode located at position x), the resultant 
dimensionless differential charge- conservation equation may be written as 


d 2 E' 

lx 2 


„E' 


-e 


[1 + e E '] 2 


( 3 ) 


where the dimensionless groupings are defined in the nomenclature. The two requisite bound- 
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Figure Is Schematic illuetration of thin-layer cell, 
(a) Rectangular (b) Circular 
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ary conditions reflect the insulating wall at the far end of the working electrode (equation 
4 a), and the pure ohmic voltage loss from the reference electrode to the leading edge of the 
working electrode (equation 4b). 


dE* 

dX 


0, X = 1 


(4) 


E* 


EIpp + 7 


dE* 

dX 


X 


0 


(4) 


Upon specification of the applied voltage, E* pp , the reference electrode location 7 , and the 
parameter 0, the nonlinear ODE was numerically integrated using the technique given by 
Newman [4]. With the potential distribution available, the cell current is found by numerical 
integration of equation 2 along the electrode surface. If the current is normalized by the 
maximum at E 0< in the absence of resistance effects, the result is defined as I and may be 
calculated from 


I 


=4 / 


a E* 


[i+>t 


d X 


(5) 


If a circular working electrode is considered (Figure lb), the dimensionless charge-conservation 
equation may be written as 


L±(x*E\ = s 

X dX V dX ) [1 + e B *] 2 


( 6 ) 


with boundary conditions 


dE * 
dX 


= 0 , 


A' - 0 


(7) 
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X = 1 


( 7 ) 


I?' = E* 

&PV 


In 7 


dl E* 
dX 


Numerical integration of equation ( 6 ) using Newman’s technique enabled a calculation of the 
integral current as a function of the applied voltage and the two parameters 9 and 7 . 


RESULTS AND DISCUSSION 

Calculations using the pore-analog model result in voltammograms identical to those 
reported by Hinman et al. from a network model. For a circular electrode, Figure 2 illustrates 
the shape of the voltammogram for the anodic portion of the scan ( 9 < 0). On Figure 2a is 
shown the dimensionless current vs applied voltage for various values of the parameter 0 with 
the reference electrode located at 7 = 1. Figure 2b presents the same calculation for \9\ = 10, 
but with the reference electrode located at various positions upstream from the leading edge 
of the working electrode. The parameter 9 is a measure of the relative importance of ohmic 
resistance. For |0| much greater than one, the voltammogram becomes distorted; but if |0| 
is 0(1) or less, the voltammogram assumes the symmetric Nernstian behavior. For example, 
even at |0| = 10 with 7 = 1, the peak potential is displaced only 8 mV from E (25 C, n - 
1 ); whereas at |0| = 100 the displacement is 43 mV, and the wave is severly distorted. As the 
reference electrode is moved farther from the leading edge (Figure 2b), the peak displacement 
and wave distortion become more pronounced at any value of 9. 

The peak current is dependent upon 0 but not 7 ; whereas the peak potential is sensitive 
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Figure 2b. 


to both. The identical behavior is found in the rectangular electrode , but the current-wave 
distortion occurs at lower |0| values in this geometry. 

A large number of calculations were performed to generate a relationship between the 
peak potential and current and the two parameters » and y. An empirical relationship was fit 
to the data as reported below 

ip., ffi (8) 

1 pl " 1 /c(7) + I*I/*(7) 

where b(-y) and 0(7) are polynomial functions of 7 as given by 
Rectangular electrode: 


6(7) = 6.9811 f 26.62997 - 0.066057 2 + 7.084 x 10 4 7 3 ( 9 ) 

c ( 7 ) = 0.09468 + 0.26217 + 1.694 x 10 -3 7 2 ~ 1314 x 10 * 7 * 

Circular electrode: 

6(7) = -23.2683 + 29.42297 - 3.3209 7 2 + 0.1342 7 3 (H) 

c (~) = -0.0765 + 0.13187 - O.OI82I7 2 + 9.183 x 10'V ( 12 ) 


Figure 3 illustrates these correlations for E; in comparison to the numerically calculated values. 
The discrepancy between the two is seen to be the largest when the reference electrode is placed 
at the leading edge and \6\ is less than 1. However, at these values of 7 and 9, wave distortion 
would be difficult to measure experimentally and is not important. 

The peak current, I p , was also empirically correlated with 9 for each electrode geometry 
and is given below. 
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Rectangular cell : 

r P 

= 1 - 6.237 x 1O- 3 !0 | 0 - 961 

(13) 

Circular cell : 

h 

= 1 - 9.724 x lO-IOI 1 - 219 

(14) 


The correlations may be used, for example, to estimate the size of an electrode which 
can be used in a particular media to maintain the wave-distortion at an acceptable level, or 
to anticipate the distortion of the wave due to reference electrode placement. It must be 
emphasized that these relationships are empirical and should be used with caution outside of 
the range used to establish the equations: 0.1 £ |0| < 10 2 ; 0 < 7 < 10. This is not a severe 
limitation considering that at low | 0 | the distortion is unmeasurable, and no experiments 
would likely be earned out at | 0 | values much larger than 100 where the current wave is 

grossly distorted. 
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NOTATION 


D 

E 


E° 


E* 


E ; 


K PP 


F 


h 


I 


Ip 

L 

n 

r 

R 


initial concentration of oxidized species of redox couple, mol/cm 

diameter of circular electrode, cm 

working electrode voltage relative to reference, V 

formal potential of redox couple, V 

dimensionless voltage, nF(E-E° )/RT 

dimensionless voltage at peak current 

dimensionless voltage applied at reference electrode 

Faraday’s constant, C/equiv 

cell gap height, cm 

reaction current density, Eqn. 2, A/cm 2 
solution current density, A/cm 2 
dimensionless current: 

rectangular cell-/. 1 [(nFfuCUhL/iRT] 

circular cell-/,'’'’ ri„<ir/ |(nP)VC;,(D/2)V8«r] 

dimensionless peak current 
rectangular electrode length, cm 

number of e~ in redox couple 

radial coordinate, cm 

gas constant 
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RflE 

radial position of reference electrode in circular cell, cm 

T 

temperature, °K 

X 

rectilinear coordinate, cm 

XflE 

position of reference electrode in rectangular cell, cm 

X 

dimensionless coordinate: 
rectangular cell — x/L 
circular cell — 2r/D 

Greek 

7 

dimensionless location of reference electrode: 
rectangular cell — x/^/L 
circular cell 2R/*jc;/D 

K 

solution conductivity, (ohm*cm) -1 

V 

negative of voltage scan rate, — dE/dt, v/s 


6 dimensionless resistance: 

rectangular cell — (^r) 
circular cell — (^r) tmFCat ^ D ^ 2 ^ 
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ABSTRACT 


A cylindrical tube is used as the basis for a two-dimensional mathematical 
model to calculate the impedance of a flooded porous electrode. The model incorporates 
charge-transfer, mass-transfer, and ohmic resistances to obtain the axial and radial de- 
pendencies of the concentration and potential profiles. A linearized Butler- Volmer kinetic 
expression for a simple redox reaction O + e^ R is used, in conjunction with analytical 
expressions for the surface concentration and overpotential, to compute the open-circuit 
impedance. The results of the two-dimensional model, which omits double-layer charg- 
ing, are compared with the results of a more standardly applied one-dimensional model, 
in which radial variations are neglected, with and without double-layer charging. The 
simpler- to- apply one-dimensional model is found to be satisfactory when mass-transfer 
and ohmic resistances are small with respect to charge-transfer resistance. The omission 
of double-layer charging does not introduce any error into the two-dimensional model in 
the frequency range in which capacitive-like effects are caused by mass-transfer limitations. 
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1 INTRODUCTION 


Impedance spectroscopy is a powerful tool for analyzing the kinetics of electrode reac- 
tions, but its application to porous electrode systems has been hindered by the lack of 
a comprehensive model for analyzing the frequency-response data. This paper presents a 
two-dimensioned mathematical model of mass and charge transport in a tube-analog model 
of a pore in a porous electrode and has been used to predict the current density distribution 
in response to both a steady and a sinusoidally oscillating potential for a redox reaction 
O + e ft. The details of the model are described below, along with sample results for 
the open-cirCuit impedance and a comparison of results from the two-dimensional model 
with a more simplified one-dimensional model. The details of the steady-state model are 
presented elsewhere [1] . 

Considerable effort has been invested in the analysis of impedance measure- 
ments for porous electrodes; however, mathematical complexity has generally limited the 
development to one-dimensional models. One of the early workers in this area was Win- 
sel [2] who solved the two-dimensional Laplace equation for the potential distribution in 
a cylindrical pore given an arbitrary expression for the current at the pore wall. He pro- 
posed a constant interfacial impedance comprised of a sum of parallel impedances arising 
from double-layer capacitance and coupled diffusion and reaction impedances. A linearized 
Butler- Volmer expression was used to describe the current-overpotential relation. Winsel 
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simplified the problem by approximating the radial diffusion of reactants and products as 
simply linear diffusion to and from a flat plane, thereby limiting the model to large pore 
diameters and high frequencies. 

Darby [3,4] developed a model for the diffusional impedance based on a radially- 
averaged concentration while ignoring ohmic drop within the pore. De Levie [5] considered 
the case in which only kinetic resistance was significant. Keddam et al. [6] proposed 
a one-dimensional model incorporating both charge- and mass-transfer resistance while 
ignoring the effects of ohmic resistance. One of the most ambitious models was that 
of Rangarajan [7] who developed a one- dimensional model for the impedance in a pore 
including charge- and mass-transfer resistance extending from within the pore to a distance 
5 outside of the pore. A linearized Butler- Volmer expression was used to describe the 
current-overpotential relation. Reiser et al. [8] incorporated a variable-pore diameter into 
a one-dimensional microscopic model of a pore. Cachet and Wiart [9] developed a one- 
dimensional model similar to that of Rangarajan, extending it to conditions far from the 

equilibrium potential. 

Most of these models are based on a one-dimensional description of the con- 
centration or potential field within the pore. The present work follows the course first 
set by Winsel [2]; i.e., development of a two-dimensional model of the concentration and 
potential fields in a cylindrical pore. In the present case, the local impedance is a function 
of distance from the pore entrance; however, for the sake of simplicity, double-layer capac- 
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itance has been ignored. (Note, however, that an analytical expression has been obtained 
to incorporate the effects of double-layer charging [see Reference [1]]). 

The details of the two-dimensional model are presented below, along with the 
description of a one-dimensional model which includes double-layer charging. A sensi- 
tivity analysis is presented to describe the results of the two-dimensional model under a 
variety of conditions. A comparison of the impedance predictions from the one- and two- 
dimensional models indicates the conditions under which the approximations in the one- 
dimensional model fail. Finally, impedance predictions from the one-dimensional model 
including double-layer charging are compared with those from the two-dimensional model 
to assess the impact of neglecting double-layer charging from the two-dimensional model. 

2 THEORETICAL DEVELOPMENT 

2.1 Assumptions 

For the purpose of model development, a porous electrode is assumed to consist of a solid 
phase permeated by a series of identical, non-interconnected cylindrical pores that are filled 
with an electrolyte solution. Since all tubes are equivalent in this idealized electrode, the 
behavior of the complete electrode will be the same as that predicted for a single tube. 
Alternatively, it is possible to predict the performance for several different tube sizes and 
then combine the results based on a given pore-size frequency distribution, as suggested 
by Winsel [2] and de Levie [5]. 
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A cylindrical tube of diameter d extends from z = 0 to a distance L within the 
electrode. At the pore aperture (z = 0), the tube is exposed to the bulk solution containing 
both oxidized species 0 and reduced species R. The potential at the inlet of the tube is 
potentiostatically controlled. The inner wall of the tube is electroactive and the end wall 
(z = L ) is inactive. Other assumptions are: 

1. The solid phase is isopotential. 

2. The reaction is a simple redox reaction with both reactant and product soluble 
(0 + e w* R). 

3. No other processes limit or alter the rates of mass transfer and charge transfer (i.e. 
no adsorption or desorption and no preceding or following reactions). 

4. The electrolyte is well supported so that the solution phase potential is governed by 
the Laplace equation and migration is insignificant. 

5. Diffusion is the only mass-transfer mechanism. 

6. There is no uncompensated resistance from the reference electrode to the pore aper- 
ture. 

7. The reactant and product species have equal diffusivities ( T>o — T)r = V). 

8. The electrolyte reservoir at the pore aperture is well mixed so that there is no mass- 
transfer resistance from the bulk solution to tube inlet. 

9. The potential and reagent concentrations are radially uniform at the mouth of the 
tube. 


Assumptions 1-5 are common to most models of porous electrodes, although some models, 
such as that of Grens and Tobias [10,11], include the effects of migration and variable 
conductivity. Grens [12] addresses the importance of these and other assumptions for the 
prediction of steady-state current density with a one-dimensional model. For the systems to 
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be considered here (i.e. simple redox reactions in well-supported electrolyte) it is reasonable 
to neglect those effects. Assumptions 6-7 limit the scope of the model slightly, however, 
they are invoked for the sake of mathematical convenience. If different diffusivities were 
included for oxidized and reduced species, then separate concentration profiles would result 
for each and the amount of “bookkeeping” required in the calculations would be increased. 
The greatest limitations in the proposed model arise from Assumptions 8 and 9. An 
approximate boundary condition could be invoked to eliminate the need for Assumption 8 
as was done recently in a similar analysis [13]. A more difficult problem arises in trying to 
eliminate Assumption 9. In reality, the concentration and potential gradients at the inlet 
of the tube will be radially dependent. The concentration gradient will be greatest near 
the wall of the tube and decrease towards the centerline of the tube. A rigorous treatment 
of this problem would require the simultaneous solution of the diffusion equations inside 
and outside of the tube — a fairly intractable problem. 

2.2 Solution to Transport Equations 

The development of the impedance model begins with the general time-dependent de- 
scription of mass and charge transfer in a cylindrical tube. Once the transport equations 
have been solved, then the specific case of impedance measurements can be simulated by 
application of a sinusoidal potential at the aperture. 

The concentration profile for both species of the redox couple is governed by 
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Fick’s second law 


( 1 ) 


dCp 

dt 


W 2 Co 


where Co is the 0 species concentration and D is the diffusivity. The concentration is 
initially uniform throughout the tube and equal to the bulk concentration. The concentra- 
tion at the inlet of the tube is constant at the bulk concentration and there is no diffusion 
through the end wall of the tube. In the radial direction, the concentration profile is sym- 
metrical about the centerline of the tube and the flux at the tube wall is proportional to 
the reaction rate. These constraints can be expressed mathematically as 


(la) t — 0; Co(r,z) — Cq 

(lb) z = 0;Co(r,0) = CS (Id) r = 0; ?§f- = 0 

(lc) z = L\ - 0 (le) r = d/2; -V&? = j(z, t ) 

where Cq is the bulk concentration of the species and j(z,t ) is the molar reaction flux at 
the wall of the tube. 

Laplace’s equation governs the solution potential $2 

0 = V 2 *2 (2) 

As with the concentration profile, the potential is symmetrical about the centerline of the 
tube and the current density at the tube wall is the sum of the reaction current density 
and double-layer charging. It is assumed that no current flows through the end wall of the 
tube and at the inlet the potential varies according to a known time-dependent function. 
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These constraints can be expressed mathematically as 


(2a) * = 0; tj = $!(() (2c) r = 0; = 0 

(2b) z=l;!£= 0 (2d) r = d/2; = i.(z, 1) + 

where $£(t) is the applied potential at the mouth of the tube, k is the solution conductivity, 
i n {z,t) is the reaction current (related to j(z^t) of Equation le by Faraday’s law), C& is 
the double-layer capacitance, and $1 is the solid phase potential (assumed constant). 

The concentration and potential profiles are obtained by first solving Equa- 
tions 1 and 2 for a delta function distribution both in space and time for j(z,t) and i n (z,f). 
A separation of variables and Laplace transform procedure was used to find these solutions 
as described in Reference [1]. The resulting dimensionless concentration and overpotential 
are indicated by 0(i?, Z y Zq^t) and fj(R } Z, Zq,t), respectively. A convolution (in time) and 
superposition (in space) integral was then applied to write 6 and rj for any arbitrary j(z , t) 
and t n (z, t) as 

0(R } Z , t) = f r 6{R > Z, i - r) J{Z oy r) dr dZ 0 (3) 

Jo Jo 

T](R,Zj)= [ l r Tj{R,Z,Z 0 ,t- t)1(Z 0 ,t) dr dZ a (4) 

Jo Jo 

where c7^(2 oj t) and T(Z 0 ,t) are the dimensionless analogs of j(z,t) and i n ( 2 ,f) as defined 
in the Nomenclature section. Note that the dimensionless time scale is different in the two 
equations; a characteristic diffusion time (tp) an d a double-layer charging time {tdi) are 
used. Further simplification is possible by recognizing that only the quantities at the tube 
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wall are of interest. Therefore, we define 


6 w {Z,t) = 8{R = l 1 Z,t) 


r u(Zj) = r,(R=l,Zj) 


and write 


= -4± f f 

n=l J0 J0 I 


X 


1 + ■£ 

J = 1 


dr dZ 0 


and 
Vu ,{Z,t) 


, ' ( ‘' ) + i?w)' in ^ z /o f, " (TK ’' <, '' r)rfr 

--f; S in p n Z f f \L{Z 0 ,T)smP n Zoe- a ' { ' t - r) ] dr dZ, 

y n=l 0 0 


( 5 ) 


( 6 ) 


where the analytical expression for f\ and 6 have been used and f3 n — (n \)^ H > are 
the zeroes of the Bessel function Ji(5), Io{u) and Ii{u) are modified Bessel functions, and 
<T n = PnIi{Pn)/Io(Pn)- No further integration of Equations 5 and 6 can be carried out until 
functional dependencies are specified for J^(Zo,t ) andX(Zo,t). At steady-state conditions, 
equations 5 and 6 are linearily related because both the potential and concentration profile 
are governed by the Laplace equation and have similar boundary conditions; hence, only 
one is needed to calculate the reaction profiles in this circumstance. 


In summary, Equations 5 and 6 give the expressions for the concentration and 
overpotential at the tube wall, respectively, and were generated by the solution to the 
two-dimensional fields subject to the assumptions and boundary conditions listed above. 


29 


Note that neither the overpotential function rj°[t ), the mass flux J(Z,t), nor the reaction 
current T(Z, t ) at the wall have been specified in these equations and can be any arbitrary 
function. Of course, there is the requirement that mass flux and reaction current must be 
proportional, according to Faraday’s law. In the most general case, the electrochemical 
reaction rate is a function of the concentration and overpotential at the tube wall through 
the kinetic expression. This introduces two levels of complication in the use of Equations 5 
and 6. First, the kinetic equations are typically non-linear and second, the equations are 
coupled, since J(Z, i ) and 1(Z , t ) are both functions of 0 W and r) w . An added complication 
arises from the two different characteristic times used to make these equations dimension- 
less; nevertheless, Equations 5 and 6 can be evaluated numerically. A considerable savings 
of computation time is realized by linearization of the kinetic expression. For the purpose 
of model development, it is further assumed that the capacitance of the double-layer is 
negligible, so that only faradaic current need be considered. This assumption simplifies 
the analysis because the double-layer charging time scale is eliminated; however, it has 
the disadvantage that it eliminates double-layer charging from the impedance model and 
consequently limits the frequency range over which the model is applicable. Additional 
analysis has been completed to include double-layer charging in the impedance model at 
the expense of added complexity [1]; however, the results presented below do not include 
it. 


Elimination of double-layer charging does not affect the expression for the con- 
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centration profile (Equation 5), however, the overpotential equation will be different. The 
modification is accomplished by dropping the second term on the right side of boundary 
condition 2d. With this simplification, the solution to Equation 2 at the pore wall is 

= rf(l) - rj(Z„,f)sin/3„Z,dZ„ (7) 

7 ^ Pnh\.Pn) 

This result is equivalent to Equation 9 of Winsel [2], in which double-layer charging is also 
ignored. With this result, the time dependence is now entirely in terms of the diffusion 
time scale. 


2.3 Impedance Calculation 

In principle, Equations 5 and 7 could be used to calculate the current- voltage-time behavior 
for any arbitrary imposed voltage and kinetic equation. However, the numerical evaluation 
of these equations would be quite difficult since it involves the double integration of the 
product of infinite sums (Equation 5). In addition, for a non-linear kinetic equation (e.g. 
the Butler- Volmer expression) an iterative procedure would have to be applied at each 
time step in the calculations. Fortunately, this calculational burden can be avoided in 
calculating the impedance at the open-circuit potential as shown below. 

For small-amplitude oscillations about an arbitrary steady-state potential, the 
dimensionless mass flux can be approximated by a first-order Taylor- series expansion 

J{0,V) = Jm + <*i(0 ~ Om) + « 2 (q - r, m ) (8) 
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where J m is the mass flux at the steady-state potential, and cti and ctj are the partial 
derivatives of J with respect to 6 and tj , respectively. For the case of oscillation about the 
open circuit potential J m — 9 m = rj m = 0; however, for the sake of generality we will retain 
the form 

= *J° + oiid w ^Zyt'j + at2»/ w (^,t ) (9) 

where incorporates all of the steady-state terms. A similar equation can be written for 
the reaction current density 

= + a 3 9 w (Z, t) + cnrj w ^Z,i) (10) 

where T° is defined in the same manner as J°. 


The definitions of the Taylor-series coefficients a,- depend on the kinetic expres- 
sion. For a system governed by the Butler- Volmer expression, the explicit expressions for 
<*i and a 2 are 






c% 


'R 


■) 


* (s). ■ - r ' 


«*(i - 0 m7 %)e~*"‘ + a c ( 1 + 0 m )e-w 


(11) 

( 12 ) 


where T 2 is a constant that relates the reaction current to mass flux and the quantities 
and e*c are the anodic and cathodic transfer coefficients, respectively. Given that the 
reaction current and mass flux are related by Faraday’s law, the ctj are related 
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For oscillations about the open-circuit potential, rj m = 0 and 6 m = 0 so the a; 
are constant. However, for conditions other than the open-circuit potential, tj*, 0 my and 
hence the a* are functions of position by Equations 11 and 12, Consequently, the quantities 
J° and J°, are also functions of position Z and the analysis that is required becomes much 
more difficult than that which is presented below. Therefore, we shall proceed on the basis 
that and the a* are constants, assuming that an appropriate (approximate) average 

values can be used when r} m ^ 0. 

Equation 10 can be integrated over the pore length to yield the total current 
once expressions are found for 0 W and rj wt The appropriate functional forms are based on 
an expansion in the eigenfunctions of the partial differential equations (Eqns. 1 and 2). In 
this case 

e w {Z,i) = f> n (t)sin/3 n Z (15) 

n=l 

T} w (Z,t) = T]°(i) + ^2 On(<) sin PnZ (16) 

n=l 

The coefficients a n (t ) and b n (t ) are unknown functions of time which can be determined by 
substituting the linearized current and flux (Eqns. 9 and 10) and the assumed expressions 
for 9 W and t} w (Eqns. 15 and 16) into Equations 5 and 7. Details of the manipulations are 
summarized in the Appendix. The total dimensionless current density is thus calculated 
as 

It = f \I° + C*30u>(^) + 04 T 7t»('^)] d-Z 

J 0 
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( 17 ) 


= t p* + )) + Ef M"(0 + «<“»(<)] 

n=l Pn 

where explicit expressions for a„(f) and ^n(^) are given by Equations 55 and 56 in the 
Appendix. Taking rj„ = 0, Equation 17 can be rearranged to yield 

Jr = ,4cosu;t + B sinwt 

= \/.4 2 + B 2 cos(u;t — <f>) 

where 

(j> = tan -1 (5/.A) (20) 


(18) 

(19) 


A = 


+ 


2a 4 ~ 
7- — 2^ 


7 /?«(«< + 

i .ffeVv-i 


’fcn 


(*k») S L + U 7 


( 21 ) 


2a 4 _ ~/p n V- 1 

b UT] C z^ hr ^ 7vT. 

7 a l „= 1 VP"/ tea 2 '** 


« fcn 


(«ta) **„ + 


(22) 




(l “ — ) 

\ a 4 + a n J 


(23) 


D-M = (^ + t)-( 


J 1 (Xfcn) 

X *" _ 2x^/ Mxk «) 


(24) 




= \/ 5 fcn + P 7 


(25) 
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^o(x) an d Ji(x) are Bessel functions and are the zeros of the function D(s) given 
by equation 52 in the Appendix* Finally, the dimensionless impedance Z can then be 
calculated as 


y = Z = |Z|e* 

= ZR e +jZl m 

where 

|z| = vCT 


(26) 

(27) 


(28) 


Several observations can be made concerning these equations. First, values of 
functions A and B, and hence Z, are determined by the three dimensionless groups 07 , a*, 
and 7 (length-to-radius ratio) and the dimensionless frequency u>, which is related to the 
dimensional frequency / by the expression 


ijj = 2ir/ 


my 

V 


(29) 


The quantities 07 and a 4 can be interpreted as 

_ mass-transfer resistance 
1 charge-transfer resistance 


_ ohmic resistance 

4 charge-transfer resistance 

These definitions will be useful in the interpretation of the results presented in the Discus- 
sion section. 
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3 ONE-DIMENSIONAL IMPEDANCE MODEL 


The two-dimensional model presented above is a rigorous treatment of transport in a 
cylindrical tube, however it is incomplete due to the omission of double-layer charging 
effects. In order to estimate the seriousness of that omission, a simpler, one-dimensional 
model was developed with double-layer charging included. A comparison of impedance 
predicted for the two different models will provide insight into the importance of double- 
layer charging. The comparison will also provide a means to assess when the radial effects 
neglected by the more commonly used one-dimensional models are important. 

The one-dimensional model material- and charge-conservation equations in a 
dimensionless format are 


96 

di 


d 2 e 
dZ 2 
d 2 e 
az 2 


-Iu 


— 2a.\6 — 2ot2Tj 


Z = 0; 6 = 0 


Z = 7! 



(30) 


and 


dt 



+ 4\on6 — 2\ot4T) 


Z = 0; t) = e jQ{ 



(31) 


where X\& is the reaction current at position Z and the concentration and overpotential 
are radially-averaged quantities. For convenience, Equations 30 and 31 were derived with 
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the same coefficients a{ as appear in the two-dimensional model. An important feature of 
Equation 31 is the inclusion of double-layer charging, which gives rise to the parameter 


K {df 2 ) 

2 CmD 


(32) 


where C < u is the double-layer capacitance (per unit area) of the interface. In effect, A is 
the ratio of the diffusion time scale to the double-layer charging time scale. 


The solution of these equations is straightforward and is found as 

0 = G a (Z)e J * F 
q = G 3 {Z)e» l 


where 


0,(Z) 


2a 2 coshfr^Z — 7)] _ cosh[r 2 (Z — 7)) 
rj — r| cosh(r a 7) cosh(r 2 7) 


(33) 

(34) 


(35) 


G 2 (Z) = 


ri - n 


{r\ - (2 a t + ju > )) 


cosh [77 (Z — 7)] 
cosh(rj7) 


/ , ^ cosh[r 2 (Z - 7)] 

-(rS-(2a 1+J uO) cosh(r2T )~ 


(36) 


and 


r, = •Jl [* + 
r, = \Jl [* - 

f = 2(o, + a.) + fl + 1) ju 


(37) 

(38) 


(39) 
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( 40 ) 


/2a, \ 

II = 4ct4(ai 4- 2a 2 ) — y + + 2a < J jv 

The average dimensionless current density is the sum of the reaction current and the 
double-layer charging current 

I ld = - 2 T l a 1 ^e{Z)dZ- 2 T 1 a 2 £rj(Z)dZ + T 1 £ \^J dZ (41) 

The impedance is computed as the ratio of the potential to the resulting current. After 
manipulation the complete expression for Z becomes 

, J" 1 (7 -\ 2a, ju>\ tanh(r l7 ) 

2«, + t J ( r .-H J 4 

, 1 i , ..\ 2 a, tanh(r j7 ) 

- U 2 ° 4 + tJ ( r = - H - a j — — 

Note that the one-dimensional current density given by Equation 41 is defined 
on the basis of the aperture area of the pore (?rd 2 /4) whereas the current density resulting 
from the two-dimensional model is based on the total internal area of the pore (ndL). 
Consequently, the dimensionless impedance computed by the two-dimensional model must 
be divided by 2 in order to compare with the results of the one-dimensional model. 

4 RESULTS AND DISCUSSION 

The approach taken in this study was to conduct a parametric evaluation of the two- 
dimensional impedance model for 27 different combinations of the three governing dimen- 
sionless groups: a,, 04, and 7. Three values of aspect ratio (7 = 5,20,50) and three values 
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of the resistance parameters (07,0:4 = 0.1,1.0,10.0) were used in the calculations. These 
values span the range from short pores (7 = 5) to fairly long pores (7 = 50), and a wide 
range of resistance ratios. A summary of the results of the parametric evaluation will be 
presented below, followed by a comparison of the results from the two-dimensional model 
with those obtained from the one- dimensional model in the limit of negligible double-layer 
capacitance. The comparison of the two models will indicate when it is appropriate to 
ignore radial dependencies of concentration and overpotential. Lastly, the results of the 
two-dimensional model will be compared with those obtained from the one-dimensional 
model when double-layer charging is included so as to assess the importance of its omis- 
sion from the two-dimensional model. Complete details on the results of the parametric 
study and further discussion of the numerical calculational procedure can be found in 
Reference [1]. 

4.1 Parametric Evaluation of the Two-Dimensional Model 

4.1.1 Aspect Ratio 

The effect of aspect ratio (i.e., pore length) is illustrated in the Nyquist plot in Figure 1. 
Dimensional impedances can be obtained by multiplying the dimensionless quantities by 
d/2K (with units of fl-cm 2 ). A total of 36 points were calculated for each curve, spanning 
12 decades of frequency. A spline routine was used to draw smooth curves between the 
points. A limited number of calculated points are shown explicitly in the Figure with the 
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Figure 1: Dimensionless impedance as a function of aspect ratio (oj = = 0.1). 
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corresponding dimensionless frequencies. The three curves in the Figure correspond to 
short, medium, and long pores. 

The model shows that the effective capacitance and resistance are both affected 
by pore length. The maximum value of Zj m decreases with increasing aspect ratio as does 
the low frequency limit of Zr € . The curves for 7 = 20 and 7 = 50 converge at a frequency 
of u> = 0.3, suggesting that both behave as semi-infinite pores above this frequency. That 
is, the impedance predicted at frequencies above u> = 0.3 for a pore with an aspect ratio 
of 20 is the same as that expected for a semi-infinite pore when ai = 04 = 0.1. 

An interesting feature of Figure 1 is the inflection that occurs in the curves for 
7 = 20 (between u> = 0.03 and 0.3) and 7 = 50 (between & = 0.003 and 0.03). For some 
cases reported in Reference [1] (e.g., when c*i = <*4 = 0.1 and 7 = 50), a local minimum 
occurs, giving the appearance of two separate loops that are added together. The feature 
is believed to be real (as opposed to a numerical aberration), since similar behavior is 
predicted by the one-dimensional model of Rangarajan [7] and the one-dimensional model 
described in Section 3. Presumably the inflection and the local minima are manifestations 
of the complex, non-spatially uniform interacions between kinetic, ohmic and mass-transfer 
effects. 
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4.1.2 Mass-Transfer Resistance 


The effect of mass-transfer resistance on the impedance (as measured by the parameter 
aj) is shown in Figure 2 for a pore aspect ratio of 20 when ohmic resistance is small 
with respect to the charge transfer resistance (a< = 0.1). Note that the smallest loop in 
Figure 2 (a! = 0.1) is identical to the loop for 7 = 20 shown in Figure 1. As expected, 
the maximum values of Zr c and Z/ m increase with increasing mass-transfer resistance 
(i.e. cti). However, the high frequency limit of Zr c is independent of the mass-transfer 
resistance. The reason for this behavior is that in the limit as or —> 00, the reaction 
is confined to a thin region near the inlet of the pore where the concentration is fixed 
through the boundary condition, hence, the high frequency limit of Zr k is not dependent 
upon mass-transfer resistance. When ax = 10, neither the low-frequency limit of Zn e nor 
the maximum value of Zj m is evident, since they occur at frequencies below the lowest 
value used in the calculations. 

4.1.3 Ohmic Resistance 

The effect of ohmic resistance on the impedance of a tubular pore is illustrated in Figure 3. 
The curves in the Figure for a< = 1 and 04 = 10 are scaled by the factors 10 and 50, 
respectively, to make them more easily visible on the plot. In fact, the Zr { values for 
04 = 10 are actually less than those for 04 = 1, which, in turn, are less than those for 
04 = 0.1. This apparently anomalous behavior is actually an artifact of the manner in 
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Figure 3: Dimensionless impedance as a function of ohmic resistance. Impedance values 
computed with a, = 1 and o< = 10 are scaled by factors of 10 and 50, respectively. (a x = 0.1, ' 
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which Zfi e is made dimensionless. Specifically, the dimensional impedance is multiplied by 
2k/ d. Therefore, as the conductivity decreases, the product of the dimensional impedance 
and conductivity also decreases. The trend of the capacitance is intuitive: as the ohmic 
resistance increases, the resulting faradaic and capacitive current decrease. 

4.2 Comparison of One- and Two-Dimensional Models for Cal- 
culation of the Impedance 

The major advantage of the two-dimensional model is that it can be applied in cases 
where the radially averaged quantities inherent in the one-dimensional models result in 
poor approximations. Furthermore, with the two-dimensional model available, it can be 
used to assess the error in one-dimensional approximations. The major disadvantage of 
the two-dimensional model is that it requires a greater amount of computational effort. 
The purpose of this section is to assess when that effort is necessary and when the one- 
dimensional approximation is sufficient. For the results described here, the double-layer 
capacitance in the one-dimensional model was set to zero. 

4.2.1 Parametric Test of One-Dimensional Impedance Model 

The one- dimensional model is expected to fail when radial gradients become large (i.e. 
when mass-transfer and ohmic resistances are large with respect to the charge-transfer 
resistance). The error in the impedance resulting from the neglect of radial gradients in 
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concentration and overpotential profiles is illustrated in Figures 4-6, which compare 
the results of the two-dimensional model with those from the one-dimensional model. The 
data presented in Figure 2 are used as the example for comparison. As expected, when 
mass-transfer resistance is small with respect to charge-transfer resistance (i.e. a i = 0.1, 
Figure 4), the one-dimensional model mimics the two-dimensional model quite well. The 
peak value of Zj m appears to be the same for the two but the real component of the 
impedance is slightly shifted for the one-dimensional model. At higher values of mass- 
transfer resistance, (Figures 5 and 6), the absolute difference between the predictions 
for the two models increases. Whereas the maximum value of Z/ m is approximately the 
same for the two curves in Figure 4, the maximum value for the one-dimensional model 
is lower as the mass-transfer resistance increases. It may be expected that an increase in 
the mass-transfer resistance would yield an increase in the capacitive component of the 
impedance. However the one-dimensional model tends to dampen the effects of mass- 
transfer resistance by equating the wall concentration of reactant species with the average 
value at that location. There is a small difference between the values Z Re {u? -► oo), however 
it is not evident at the scale shown in the Figures. Although the curves in Figure 6 do not 
show the low frequency loop, it is evident that the same trends will apply. 

The effect of ohmic resistance on the one-dimensional model is illustrated in 
Figure 7, which shows a log-log plot of the impedance modulus as a function of frequency. 
The three pairs of curves in the figure are the results from the two-dimensional model 


46 




m^ 4: ( 0 ?=r=“l° f 7 ‘ = M) diCted imped “ ce fOU “ d fr ° m the 0ne - “ d ‘^dimensional 





F i g ure 5: Comparison of the predicted impedance found from the one- and two-dimensional 
models, (a, = 1 , a< = 0 . 1 , 7 = 20 ). 


48 




Figure 6: Comparison of the predicted impedance found from the one- and 
models. ( 0l = 10, a 4 = 0.1, 7 = 20) 


two-dimensional • 
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Figure 7: The modulus of the impedance as a function of frequency for 1-d (solid line) and 
2-d (dashed line) models, (oj =0.1, 7 = 20). 
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presented in Figure 3 and the corresponding results from the one-dimensional model. As 
a 4 increases, the modulus decreases and the slope of the modulus-frequency curve goes 
to zero. The zero slope is a result of the fact that the effective capacitive component of 
the impedance goes to zero. (The modulus decreasing with increasing a 4 is an artifact of 
the way in which the impedance is made dimensionless, as discussed in conjunction with 
Figure 3.) The results shown in the Figure suggest that the one-dimensional model would 
be a poor choice for analyzing data collected from a system with high ohmic resistance. 

In many applications, the dimensionless parameter a 4 will be less than 1 and 
typically on the order of 0.1 or less. Note, however, that it would be an oversight to 
ignore the effects of ohmic resistance completely. As an illustration, consider the set of 
curves shown in Figure 8, which were obtained with the one- dimensional model. At a 
value of ct 4 — 0.1, ohmic resistance is low with respect to the charge-transfer resistance. 
However, decreasing that value by half (i.e. a 4 = 0.05) produces a significant shift in both 
the real and imaginary components of the impedance. A ten-fold decrease in a 4 (from 
0.1 to 0.01) increases the maximum value Zj m roughly three-fold. In the limit of very low 
ohmic resistance, the impedance plot appears quite different from that obtained with what 
appears to be a “low” value (i.e. a 4 = 0.1). Therefore, one should use care in applying 
models that ignore the ohmic potential drop completely such as that of Keddam tt a/.[6]. 

Since both the one- and two-dimensional models account for axial gradients 
in a finite pore, it is intuitive that the agreement between the two models would-not be 
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affected by the aspect ratio and is not discussed further here. 

4.2.2 Applicable Parameter Range of the One-Dimensional Model 

The results presented above show that the assumption of radially-averaged concentrations 
and potentials can be a serious limitation to the one-dimensional model. The problem 
arises when mass-transfer resistance is high, and concentration gradients extend far from 
the wall of the tube. Likewise for the relation between ohmic resistance and the radial 
gradient of the potential. Consequently, a radially-averaged quantity can be significantly 
different from the actual value at the wall, with the result that the current predicted from 
averaged quantities will differ from that predicted from the two-dimensional model. 

The problem of radial averaging will be most significant at very low frequen- 
cies and/or when diffusion and ohmic resistances are large with respect to charge-transfer 
resistance. Based on the results of the parameter study it appears that, over the frequency 
range where double-layer charging is not important, the one-dimensional model is an ade- 
quate approximation to the two-dimensional model when ohmic resistance is low (a* < 0.1) 
and mass-transfer resistance is low-to-moderate with respect to charge- transfer resistance 

(<*i < 1). 
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4.3 Assessment of the Double-Layer Effects 


The most significant limitation of the two-dimensional model for the impedance calcula- 
tion is the omission of double-layer capacitance; however, the one-dimensional model does 
incorporate it and a comparison of the two will indicate when the omission adversely af- 
fects the results of the two-dimensional model. As described in the preceding sections, 
the one-dimensional model is not without limitations. Nevertheless, for cases where mass- 
transfer and ohmic resistances are small relative to the charge-transfer resistance (e.g. 
when ax = 04 = 0.1), the one-dimensional model is an adequate approximation to the two- 
dimensional model. Therefore, conclusions based on the behavior of the one-dimensional 
model in these parameter ranges will be applicable to the two-dimensional model. 

The effect of double-layer charging as predicted by the one-dimensional model 
is illustrated in Figure 9. The dashed line in the Figure shows the results computed 
with the two-dimensional model and the solid line was obtained with the one-dimensional 
model, including double-layer effects. For the results shown here, the ratio of diffusion to 
double-layer charging characteristic times is A = 2.27 x 10®. The two semi-circular curves 
obtained in the frequency range u> = 3 x 10 -4 to 30 are identical to the diffusion loops 
shown in Figure 4. The large loop at higher frequencies (w > 3 x 10 4 ) is a result of the 
interaction of double-layer charging with charge- transfer kinetics. 

The significance of these results can be seen more easily in Figure 10, which 
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Figure 9: A comparison of the impedance calculated from the two-dimensional model (dashed 
line) with the one-dimensional model including double-layer charging (solid line): (a, = a. = 
0.1, 7 = 20, A = 2.27 x 10 6 ). JK 








shows the modulus of the impedance (|Z|) and phase lag ( <f > ) as a function of frequency 
for the one-dimensional model with and without double-layer charging as well as the two- 
dimensional model. There are two points to be made from these results. First, in the 
low frequency range (u> < 10 4 ), the one-dimensional model with double-layer charging is 
identical to the one-dimensional model without it. The overlap of the results from the 
one-dimensional models is quite evident. The agreement between the two models in the 
low frequency range leads to the conclusion that the omission of double-layer charging does 
not introduce any error into the prediction of the diffusion loop in the impedance diagram. 
By analogy, we can conclude that the results presented for the two-dimensional model are 
also correct as presented. 

The second important result shown in the Figure is that it clearly illustrates the 
wide frequency separation of the diffusion and double-layer loops that appear in Figure 9. 
The diffusion effect is confined entirely in the frequency range below uj = 100. That is, 
at frequencies above a) = 100 mass-transfer resistance becomes insignificant in the vicinity 
of the pore wall. Conversely, the double-layer charging effect is not significant below 
frequencies of u> = 10 4 . Thus, for the case considered here, diffusion and double-layer 
charging are independent effects which can be considered separately. The fact that, the two 
curves do not interact with each other arises from the large difference in the characteristic 
time scales of the phenomena involved. Diffusion, a relatively slow process, is important 
at low frequencies and insignificant at high frequencies. Conversely, double-layer charging, 


57 


which is a relatively fast process, can only be detected at high frequencies. 

The frequency separation between the two loops in the impedance diagram is 
solely a function of the parameter A and is not affected by the parameters and 04. 
Changing the values of ai and <*4 have the effect of shifting the curves in the Bode plot 
horizontally and vertically, but the shape of the curve is not affected. Consequently, the 
parameter A can be used to judge whether or not the omission of double-layer charging 
introduces errors into the impedances predicted by the two-dimensional model. Results 
of a numerical study show that the frequency separation between the two loops will be 
at least one order of magnitude as long as A > 3000. That is, given A > 3000, then 
the double-layer charging loop will not become evident until a frequency that is 10 times 
greater than the frequency at which the diffusion loop disappears. 

5 Summary 

A tube-analog model of a pore was used to develop the two-dimensional, time-dependent 
concentration and potential profiles in the pore. Impedance results were computed by 
application of appropriate boundary conditions to the solution obtained for the general 
problem. A comparison of the two-dimensional model with a one-dimensional model based 
on the same assumptions shows that the one-dimensional model is a poor approximation 
when mass-transfer and/or ohmic resistance are dominant. In addition, it was found that 
omitting ohmic resistance altogether can introduce errors in the computed impedance. 


58 


Although the model presented above omits double-layer charging, no significant 


errors are introduced into the mass-transfer capacitance loop. 
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6 APPENDIX 


The coefficients a n (t ) and b n (t ) can be determined by substituting the linearized current 
and flux (Eqns 9 and 10) and the assumed expressions for 6 W and rf w (Eqns 15 and 16) 
into Equations 5 and 7. When the resulting equations are multiplied by sin 0 n Z o and the 
orthogonality property of the sinusoidal eigenfunctions is applied, an integral equation is 
obtained for b n (t), and a n (f) is algebraically related to 6„(f ). 


MO = -- /' (J* + + 5 (°.MO + 

7 Pn ^ 


xe-«< r - T > x 
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dr 
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a n (t ) = 
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7( a 4 + Cn) 


7T + T"M* ) 4- —T] ( t ) + ——a n (t ) 
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where <7 n = /? n /i(/? n )// 0 (/? n ). The solutions to Equations 43 and 44 are obtainable by 
Laplace transform. In Laplace domain Equations 43 and 44 are expressed as 

5„(a) = “ + «2*r) + | (aJn + ajOn) 

* (•+« + S *+«+«;) (45) 


-2 


On = 


7 ((r„ + a*) [s/3, 


J " +^6„ W +“ 4rW 


ft, 
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Equations 45 and 46 can be manipulated to yield explicit expressions for 5 n 
based only on the Laplace variable a and the physical parameters of the system; how- 
ever, reinversion of that expression would yield an integral that could only be evaluated 
numerically. A more useful approach is to invoke the identity relation 


1 00 1 
i + E 


1 l p(g) 


(47) 


x 2 fr{x 2 + S^ 2xli(x) 

This identity was derived by comparing the concentration profile obtained from the steady- 

state formulation of the governing equations with that obtained by considering the long- 
time response of Equation 5 to a step-change in the mass flux at the wall of the tube. The 
two expressions must be identical which results in equation 47. The identity relationship 
has been verified numerically [1]. After substitution of Equation 47 into Equation 45 and 
rearrangement, the resulting expression for 6 n (s) is 


6 n(s) 
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Substitution of Equation 46 into Equation 48 and rearranging leads to 


K{s) = 


D{s) 


( 49 ) 
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Z?(s) = a y» + «/,(^ Tfii) + «i (l - — ^ 2l_) I*{y/s + Jl )] (52) 


Equation 49 can be inverted according to the formula 


MO = £ 




*(•) 


z>( s ) 




fc=l 


where D'(s) is the derivative of the function D with respect to a 


®'M = 


+ 


v77ffi + -4— «, 


A(\A + /?n) 


(53) 


(54) 


and a* n are the roots D(skn) = 0. The inverted expression for b n (t ) is: 

where the convolution integral has been applied to the second term in the numerator of 
Equation 49. An important point to note about this expression is that it is independent 
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of the form of the imposed potential function The only restriction is that the 

perturbation must be small so that the assumption of linearity in the kinetic equation 
applies. 


When a sinusoidal potential function with amplitude V°c is inserted into Equa- 
tion 54, integrated, and transients are eliminated by taking the limit as f — * oo, then 




+ 


Vc 


h (T\ — V ^ 2 ^ Skn ) 

n{ ) ~ k [ S kn ' 


(a) sin Cot — Skn cos Cot) 


(55) 


The long-time limit of o n (t ) is obtained by substituting the definition of b n (t ) into Equa- 
tion 44 and rearranging: 

2 


' a n {i ) = ~ 


j(a 4 + <r n ) 


1 ° . <*4 037 

K + K”- ~ 


N 2 {s kn ) 1 

2 k D’(s kn ) s kn 


, . _j ( <*2lUT} 0 c ^ N 2 {s k n ) 

+ smwi ( — 2^ 


V 2 kl D '{*kn) U> 2 + 3 2 kn 


) 


4* cos cot 



037 « > Y' jy a( 3 fcn) Jfen V 
2 ^ kl D '( 9 kn )^ 2 + ^In) . 


(56) 


In summary, Equations 55 and 56 can be substituted into Equations 15 and 16 
to obtain the stationary-state concentration and overpotential profiles for a reaction oc- 
curring in response to a small amplitude, sinusoidally oscillating potential at the inlet of 
a cylindrical tube. 


It can be shown [1] that all s* n are reed, negative, and infinite in number 
and were found by using standard root-finding methods [14] for low-order roots and an 
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asymptotic expansion for larger roots. Thus, computing the impedance consists of finding 
the roots of the function D(skn) defined in Equation 52 and then calculating the infinite 
sums for use in finding values of A and B (Equations 21 and 22). Once these values are 
known, Zn e and Z/ m can be computed directly. 
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Nomenclature 


c dt 

Ci 

c r 

Co 

C R 

d 

V 

f 

F 

1 

Zu 

i n {z,t) 

I T 

Z„ 

j 

J 

J., 

j(z, 0 

L 

r 

l R 

I R 

t 

T 

i 

t 

I'D 

tdl 

z 

z 

z 

|Z| 

Zj m 

%Rt 


double-layer capacitance (F/cm 2 ) 

concentration of species i (mole/cm 3 ) 

concentration of species i in bulk solution (mole/cm 3 ) 

concentration of reactant 0 (mole/cm 3 ) 

concentration of product R (mole/cm 3 ) 

pore diameter (cm) 

diffusivity (cm 2 /s) 

frequency (Hz) 

Faraday’s constant 

dimensionless reaction current at the pore wall; n( z >0 

dimensionless reaction current for 1-d model, eqn. 41 
current density at the pore wall (A/cm 2 ) 
exchange current density (A/cm 2 ) 
dimensionless total current, eqn. 17 
dimensionless steady-state current 

y/=l 

dimensionless mass flux at the pore wall; 

dimensionless steady-state mass flux 
mass flux at the pore wall (mole/cm 2 s) 
length of pore (cm) 
distance in radial direction (cm) 
dimensionless radial distance; r/{d/ 2) 
ideal gas constant (J/(mole-K)) 
time (s) 

temperature (K) 

dimensionless time (on diffusion time scale); if to 
dimensionless time (on capacitance time scale); t/tji 
characteristic diffusion time; (d/2) 2 /X> 
characteristic double-layer charging time; Cji{d/2)/K 
distance in axial direction (cm) 
dimensionless axial distance; z/(d/2) 
dimensionless impedance eqn. 26 
impedance modulus, eqn. 28 

Dimensionless imaginary component of impedance, eqn. 27 
Dimensionless real component of impedance, eqn. 27 
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GREEK 


oca anodic transfer coefficient 

etc cathodic transfer coefficient 

{3 n eigenvalues of PDE, (n — 5)^/7 

7 pore aspect ratio; L/(d/2) 

T 
T 

5j zeroes of the Bessel function Ji(ti) 

77 dimensionless overpotential; «■?(*. - *0 

rj° amplitude of imposed overpotential function 

rj 9B steady-state overpotential 

rj w dimensionless overpotential at the pore wall, eqn. 7 

0 dimensionless concentration; 

d 99 steady-state concentration 

6 W dimensionless concentration at the pore wall, eqn. 5 

k solution conductivity (n _ 1 cm _1 ) 

A ratio of the diffusion time scale to the double-layer charging time scale, eqn. 32 

<7 n 0nWn)lI*{0n) 

solid potential (V) 

$ 2 solution potential (V) 

us angular frequency of imposed overpotential function (s -1 ) 

d> dimensionless frequency, eqn. 29 
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ABSTRACT 


Reported are the theoretically predicted peak potentials and currents for the cathodic 
scan of a reversible redox couple over a range of the initial concentration ratio of oxidized- 
to-reduced species from zero to infinity. The asymptotic limit for the voltammogram as 
the initial concentration ratio approaches zero, along with the infinity limit given by clas- 
sical theory, is used to develop empirical relationships which correlate the peak potential 
and current to the initial concentration ratio. Dimensionless peak values measured exper- 
imentally for different initial concentration ratios of the ferro/ferri cyanide redox couple 
are found to agree with those theoretically predicted. 
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Introduction 


Classical theoretical formulas are available to predict the peak current and potential 
at a planar electrode during linear-sweep voltammetry (LSV) for the reversible (Nerns- 
tian) reaction 

0 + ne~ *=* R (1) 

where O and R are soluble oxidized and reduced species, respectively, and n is the number 
of electrons transferred [1,2]. Matsuda and Ayabe [1] and Nicholson and Shain [2] showed 
that if the reduced species is present initially in negligible amounts and the voltage is 
swept in the negative direction, the cathodic peak current is proportional to the absolute 
value of the sweep rate to the one-half power and the corresponding peak potential is 
sweep-rate independent. Matsuda and Ayabe [1] also demonstrated that if the sweep 
begins at the equilibrium potential, E eq) the dimensionless peak current and potential 
depend only upon the product of two dimensionless ratios, 7 6 *: the ratio of the square 
root of the diffusion coefficients of O and R } 7 = \J Do / , and the initial concentration 
ratio of O and R in the bulk electrolyte, 0* = (Cq/Cr). Matsuda and Ayabe [l] presented 
calculations for the case where 7 0* — > 00. Left numerically unexamined was the influence 
of R when its concentration is not negligibly small. Farsang et al. [3], using the same 
formula as derived by Matsuda and Ayabe [1], examined the influence of 7 0* on the peak 
values over the narrow range of 0.22 < 7 Q* < 5.5 and found a shift in the potential to 
more negative values and an increase in the dimensionless peak current as 7 6 * decreased. 
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Nicholson and Sh&in [2] considered the case of starting the sweep at a non-equilibrium 
potential, 2?*, in the limit of the reduced species concentration equal to zero (0 # = oo). 
They obtained a current-potential relationship which is a function of an analogous dimen- 
sionless grouping, 7#, where 7 is the same as defined above, and 8 is the concentration 
ratio at the electrode surface given by the Nernst equation immediately after Ei is ap- 
plied. This non-equilibrium initial condition gives rise to a Cottrell current in the early 
portion of the sweep which they showed does not significantly affect the peak current 
or potential as long as 7 6 > 650. This corresponds to starting the sweep at a potential 
which is at least 200/n mV [T = 25°C) positive of the peak potential. 

In what follows, the Matsuda and Ayabe [l] equation is applied to calculate the 
peak potential and current over the range of 7#* from zero to infinity. Also presented is 
an explicit heretofore unrecognized formula for the peak values in the limit of 7$* -► 0. 
The peak values at the two limiting behaviors are used to deduce and develop an easy- 
to-apply empirical relationship which correlates the results as a function of jO* over the 
entire range. The calculations are confirmed by experimental measurements. The effect 
of mistakingly applying the Nicholson and Shain [2] formula for regions of 7 9* below its 
range of validity is also quantified. 

Peak Current and Potential Calculations 

Matsuda and Ayabe [l] solved the one-dimensional transient diffusion equation for 
both soluble redox species, and obtained the time-dependent surface concentrations for O 
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and R in terms of the reaction current. Starting the sweep at the equilibrium potential, 


6* is related to E e, by 


r = exp [(jgr) (E„ - £"')] 


( 2 ) 


where E°' is the formed potential. For the case of reversible kinetics, substituting equa- 
tion 2, the surface concentrations and the sweep rate, v, into the Nernst equation results 


in the integral equation for the current 


L 


al 

o \/ai - z 1 + 7 0 * exp(-af) 


sL(at) 


(3) 


where t is time, a = (bt) v which is positive for a cathodic sweep, and x( a *) > s the 
dimensionless current defined as 


-Ax(at) = 


i(at) 


nFAC'oiVoa) 1 /* 


(4) 


where F is Faraday’s constant, A is the electrode area, and i(at ) is the time-dependent 
reaction current which is positive for a cathodic reaction. The product (at) is the dis- 
placement of the dimensionless potential from the starting potential. 

Equation 3 differs from the analogous formula given by Nicholson and Shain (equa- 
tion 22 of [2]) in that by their assumption of Cr is equal to zero, the numerator of Z.(at) 
is equal to one, and 9* is replaced by 9 where 9 is given by 


e = exp [($) (£, - £"')] 


(5) 
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Equation 3 is an Abel integral equation [5] for *(a<) with the solution 


v/?xoo = fl + r iz 

4y / 7r /o y/at — z cosh 2 ( — 7 ** 1 ) 


( 6 ) 


The analogous current-potential relationship derived by Nicholson and Shain (equation 

33 of [2]) from the solution of equation 3 with a unity numerator for L(at ) and 9 replaced 
by 9 * is 


v^rx(at) = 


1 r at 
\/ir Jo 


dz 


V5W(1 + 7*) v /^37cOsh S 


( 7 ) 


(For the calculations presented here, the integrals in equations 6 and 7 were numeri- 
cally evaluated using the IMSL [6] quadrature subroutine DQDAGS which accommodates 
integrands having end-point singularities.) 


Although equation 6 as 7 9* — + 00 is identical to the j9 — » 00 limit of equation 7, 9* 
and 9 are not interchangeable. The leading term in equation 7 is a result of assuming 
the initial concentration of the reduced species is zero, and predicts a Cottrell current 
upon application of E { at the start of the sweep. This term can be significant when 
7 9 is not lar S e> To illustrate this point, the effect on the voltammogram of erroneously 
replacing 9 with 9'‘ in equation 7 is shown in Figure 1 with the potentials plotted relative 
to the half-wave potential, E 1/2 = E° + (Jf)ln 7 [2]. At 7 0* = 20, the dimensionless 
peak current and potential predicted from equation 6 both differ from those found using 

equation 7 (with 9 set equal to 9') by 2%, and at 7** = 2, the percentage difference is 
6% and 33%, respectively. 
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Figure 1: Predicted linear-sweep voltammograms from a planar electrode. The solid 
lines are the voltammograms generated from equation 6. The dashed lines are voltam- 
mograms generated from equation 7 with 6 replaced by 0 *. The potentials, in mV, are 
for T = 298 K. 
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The dimensionless peak current and potential may be determined as a function of 
jO* by calculating the voltammogram using equation 6 to locate graphically the peak 
position. However, a more convenient analytical procedure was derived by differentiating 
equation 6 with respect to (at) and setting the result equal to zero. The potential 
which satisfies the resulting expression is the peak potential which, when substituted 
into equation 6, specifies the peak current. Newton’s method was used to solve the 
resulting non-linear equation for the peak potential generated from this procedure. The 
same technique was applied to equation 7 (with 9 replaced by 6*) and the resulting 
dimensionless peak currents and potentials from the two predictions are compared in 
Figure 2 as a function of 7 6 *. The dimensionless peak potential, E *, is defined as 

E; = [(at), -In 7**] (8) 

where ( at) p is the dimensionless peak potential relative to the equilibrium potential. 

Nicholson and Shain [2] report their results at 7 0 = 650, and Matsuda and Ayabe 
report theirs at 7 8 * = 165, 000, and from Figure 2 we see that 6 and 0 * are interchange- 
able at these values; that is, no difference is found in the predicted peak values regardless 
whether equation 6 or 7 is used. Replacing 8* with 6 at smaller values of 6 * , however, will 
lead to differing results. Further, for values of 7 8 * < 1.60, equation 7 does not predict a 
peak current because the Cottrell current dominates the transient response. 

The dimensionless peak values (a/^Xp = 0.4463 and E* = 1.109) given by Matsuda 
and Ayabe [1] and Nicholson and Shain [2] are the limiting values as -yd* — ► 00. Limiting 
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log 70 * 

Figure 2: Dimensionless pe&k currents and potentials for linear-sweep voltammetry on 
a planar electrode as a function of 76*. The solid lines result from the use of equa- 
tion 6, and the dashed lineB result from the use of equation 7 with 6 replaced by 6*. 
For 7® ^ 1.60, equation 7 cannot predict peak values. The symbols are the peak val- 
ues which were obtained experimentally using the Fe(CN)j 4 /Fe(CN)a S redox couple. 
□ — present results; A — Tomcsanyi et al. [4] data. 
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values for the dimensionless peak current and potential can also be obtained for 70 * — » 0 
by taking the limit of equation 6 as 7 6 * approaches zero which, after some mathematical 


manipulation, results in 


dz (9) 

Following the analogous calculus procedure applied to equation 6 to calculate the peak 
potential, Newton’s method was applied to determine numerically (at),, = 0.855 in the 
limit of - yd * — ► 0. The corresponding dimensionless peak current found from equation 9 


is 0.6103. 


Empirical Correlations 

The following one parameter empirical relationship, motivated by the asymptotic 
trends seen in Figure 2, was fit to the theoretically predicted dimensionless peak currents 
from equation 6 

v^x,b n = v^Xp(O) - [v^Xp(O) - v^x,(°°)] < 10 > 

where and V^XpC 00 ) are the known dimensionless peak currents at the limits 

■yO* — ► 0 and 7 6 * — ► 00 , respectively, and m is an empirical parameter which was found 
by least-squares fit of equation 10 to the data in Figure 2. The resulting correlation is 

vfrfabn = °.6 10 3 - 0.164 [ l og^'^. ] < n ) 

The following one parameter empirical relationship, also motivated by the asymp- 
totic trends seen in Figure 2, was fit to the theoretically predicted peak potentials from 
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equation 6 


E;(-,tr) = in{[exp£;(o)r + i«p£;(oo)r} ,/m (12) 

where £*(0) and £*(00) are the known dimensionless peak potentials at the limits 7 0* —* 
0 and 70* — » 00, respectively. The parameter m was found by least-squares fit to the 
data in Figure 2. The resulting correlation is 

W) = 1- [(^?) + 3.74 (13) 

For 0 < 70* < 00, equations 11 and 13 deviate from the exact values by at most 0.1% 
and 1%, respectively, differences which are visually indistinguishable in Figure 2. The 
maximum 1% error resulting from equation 13 (which occurs in the vicinity of 70* = 3) 
translates to less than a 0.5 mV (n = 1 and T = 298 K) error which is experimentally 
insignificant. 

Experimental Results 

Linear-sweep voltammetry was conducted using the ferro/ferricyanide redox couple 
in 1.0 M KC1 at five different concentration ratios. The concentration of the ferricyanide 
was maintained at 0.010 M, and five different ferrocyanide concentrations were used: 
0.0010 M, 0.0050 M, 0.010 M, 0.025 M and 0.10 M. Cottrell experiments were run in 
order to determine 7 by taking the ratio of the product of (i-<V 2 ) for oxidation-to- 
reduction experiments, and 7 = 1.09 was obtained which agrees within 1% of that 
calculated using the diffusion coefficients for O and R given by Adams [7]. Therefore, 
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the 7#* values for the five concentrations ratios are 10.9, 2.18, 1.09, 0.436, and 0.109. 
All experiments were performed in a single- compartment cell at room temperature with 
a 1.6 mm diameter Pt-disk as the working electrode which was not perfectly flush with 
the surface, and Ag/AgCl as the reference electrode. Voltammograms were generated 
at nine different sweep rates ranging from 4 to 150 mV/6 and were run in duplicate. 
The dimensionless peak current, was obtained by plotting the dimensional peak 

current, i p , versus the absolute value of the sweep rate to the one-half power, M 1 / 2 . The 
slope, , resulting from linear regression of these data was used to calculate v^Xp 
since from equation 4 


V^X P ( 70*) 


(RT ) 1 / 2 

{nF) 3 l 2 C' 0 AV)j 2 <*M 1/2 


(14) 


1 /) 

The unknown product ADq was determined by performing a series of linear sweeps 
at which only the oxidized species was present (7#* = 00), and solving equation 14 for 
AV]! 2 using y^Xp = 0.4463. The product AVq 2 calculated in this manner was then 
used in equation 14 at the lower values of 7 0* to obtain V^Xp( 7^*)- The peak potential, 
(at) p , which is independent of the sweep rate, was read directly from the voltammogram 
for each value of 7$*, and equation 8 was used to calculate E*. 


The experimental dimensionless peak currents and potentials (open squares) are 
compared to the theoretical values in Figure 2. The experimentally determined peak 
currents deviate from those theoretically predicted by at most 1%. The measured peak 
potentials also closely agree with those predicted, except for 9 * = 10 at which a 10 
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mV deviation occurs. Also included in Figure 2 are the experimental data reported by 
Tomcsanyi et al. [4] (open triangles) from LSV run on a Pt electrode in 1.0 M KC1 
using the ferro/ferricyanide redox couple. The scatter is larger in these data, but the 
trend in ls evident. The peak potentials obtained by Tomcsanyi et al. [4] 

were reported only to the nearest 10 mV and are 30 to 40 mV higher than the theory 
predicts. 
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Abstract 


A pore in a porous matrix is modeled as a cylindrical electrode and the mass and 
charge conservation equations are solved in the context of this geometry for the simple 
redox reaction O + ne~ R when the potential at the pore aperture is ramped linearly 
in time. The resulting linear-sweep voltammograms depend on the combined effect of 
kinetic, ohmic and mass-transfer resistances, and determining the transfer coefficient, 
the standard rate constant or the number of electrons transferred is accomplished by 
correlating the peak currents and potentials to the governing system parameters. The 
correlations are established by introducing simplifications at various parameter limits. 
The first simplification is to neglect axial diffusion which allows the diffusion equation to 
be analytically solved giving the surface concentration of the redox species in terms of 
the reaction current and the solution potential at the electrode surface. In the absence 
of solution resistance, correlations are established in the limit of reversible (Nernstian) 
and irreversible kinetics which relate the peak currents and potentials to the sweep rate 
by making use of the known peak values obtained from a planar electrode and a thin- 
layer cell as the high and low sweep-rate asymptotes, respectively. To quantify ohmic 
effects in the absence of axial diffusion, an analytical solution to the Laplace equation 
for the potential is coupled to the analytical solution to the diffusion equation. The 
peak currents and potentials obtained when axial diffusion is negligible are compared to 
those obtained when the two- dimensional conservation equations are solved numerically 
in order to determine under what conditions axial diffusion affects the voltammogram. 
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1 Introduction 


Electroanalytic techniques used to study the kinetics on planar electrodes do not always 
lend themselves to studying flooded porous electrodes. For example, hydrodynamic 
methods (e.g. rotating-disk electrodes [1]) are not useful since the reaction surface is 
inside the porous matrix, and therefore controlling forced convection only affects the 
material transport to the pore mouth but not within the matrix. In contrast, linear- 
sweep voltammetry (LSV) in which the electrolyte is stagnant can be used to perform 
in situ kinetic studies on flooded porous electrodes. However, the presently available 
mathematical methodology that has been developed to determine kinetic parameters 
using LSV on a planar electrode [2, 3, 4, 5] and in a thin-layer cell [6,7,8] is not applicable 
for a pore except at limiting sweep rates. At high sweep rates the diffusion layer i6 
small relative to the pore diameter and the voltammogram has the characteristics of 
LSV on a planar electrode when ohmic resistance is negligible. For deep pores and low 
sweep rates, radial and axial concentration gradients are small for the major fraction of 
the pore volume and consequently the pore behaves like a thin-layer cell. However, at 
moderate sweep rates in the presences of ohmic resistance and axial diffusion, LSV can 
not be characterized by these two limiting cases, and the result of LSV depends on the 
combined effects of ohmic, mass-transfer and kinetic resistances. 

This paper considers two-dimensional mass and charge transport in a tube-analog 
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model of a pore in a porous matrix for the reaction 


O + ne ** R (1) 

where O and R are soluble oxidized and reduced species, respectively, and n is the 
number of electrons transferred. If the rate of charge transfer is governed by the Butler- 
Volmer kin etic expression, the peak current and potential is Bhown to be a function of six 
dimensionless parameters. Solving the coupled material and charge conservation partial 
differential equations numerically is computer intensive (on the order of an hour of CPU 
per simulation on an IBM/3090), and so carrying out the simulations to correlate the 
peak currents and potentials to the governing system parameters would be formidable. 
Therefore, the approach taken in this study is to conduct a parametric evaluation of the 
two-dimensional model by introducing simplifications at various parameter limits, and 
use a finite-element numerical solution to bridge the gap between the limiting cases. 

The first simplification is to neglect axial diffusion which allows the diffusion equa- 
tion to be analytically solved giving the surface concentration of the redox species in 
terms of the reaction current and the solution potential at the electrode surface. Fur- 
ther ass umin g solution resist ance is negligible reduces the computation time over three 
orders-of-magnitude from the finite-element solution procedure which facilitates corre- 
lating the peak currents and potentials to the governing system parameters. In the limit 
of reversible (Nemstian) and irreversible kinetics, correlations are established which re- 
late the peak currents and potential to the sweep rate by making use of the peak values 


86 



obtained from a planar electrode [2,5] and a thin-layer cell [6,7]. Analogous correlations 
were obtained by Aoki et al. [9] for LSV in a rectangular cell for reversible kinetics and 
negligible ohmic resistance, and later confirmed experimentally [10]. 

The correlations presented here, which assume solution resistance is negligible, are 
shown to be useful also when solution resistance is appreciable by a priori accounting 
for ohmic effects. The effect of solution resistance on the voltammogram is quantified 
by coupling the analytical solution to the one-dimensional, transient diffusion equation 
to the analytical expression for the solution potential at the electrode surface obtained 
by solving the two-dimensional Laplace equation [11,12]. In parameter regions in which 
numerical difficulties prevent voltammograms from being generated from these coupled 
non-linear equations, the finite-element numerical solution is used. 

These correlations can be used to determine whether the electron transfer is re- 
versible or irreversible; if the kinetics are reversible, the correlations may be used to 
estimate the number of electrons being transferred, or if the kinetics are irreversible, the 
transfer coefficient and the standard rate constant may be extract from LSV measure- 
ments. The results are not valid if material axially diffusing through the pore mouth 
contributes appreciably to the current. The influence of axial diffusion will depend on 
the sweep rate for a given pore length since the slower the sweep, the larger is the frac- 
tion of the pore volume that will be effected by axial diffusion. Again, the finite-element 
numerical solution is used to determine the parameter regions in which the correlations 
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axe not valid due to axial diffusion. 


2 Mathematical Model 

2.1 Assumptions 

For the purpose of model development, a porous electrode is assumed to be a collection 
of identical, non-interconnected cylindrical pores that are filled with electrolyte. Since 
nil pores are equivalent in this idealized electrode, the behavior of the complete electrode 
will be the same as that predicted for a single pore. Alternatively, it is possible to predict 
the performance for several different pore sizes and then combine the results based on a 
given pore-size distribution, as suggested by Winsel [11] and de Levie [13]. 

A cylindrical tube of diameter d extends from z = 0 to a distance L within the 
electrode. At the pore aperture (z = 0), the tube is exposed to the bulk solution 
containing both oxidized species O and reduced species R. The potential at the inlet of 
the tube is potentiostatically controlled. The inner wall of the tube is electroactive and 
the closed end ( z = L) is inactive. Other assumptions are: 

1. The solid phase is isopotential. 

2. The reaction is a simple redox process between soluble reacting species. 

3. No other processes limit or alter the rates of mass transfer and charge transfer (i.e. 
no adsorption or desorption and no preceding or following reactions). 

4. The electrolyte is well supported so that the solution-phase potential is governed 
by the Laplace equation. 
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5. Diffusion is the dominant mass-transfer mechanism. 

6. The diffusion coefficients for the two reacting species are equal {Vo = Vr = V). 

7. Non-faradaic current is negligible. 

8. There is no uncompensated ohmic resistance from the reference electrode to the 
pore aperture. 

9. The electrolyte in the reservoir is stagnant and the axial concentration gradient of 
the reacting species at the pore opening is proportional to the difference between 
the bulk and inlet concentrations. 

10. The potential is radially uniform at the pore aperture. 


Assumptions 1-5 are common to most models of porouB electrodes, although some, 
such as that of Grens and Tobias [14,15], include the effects of migration and variable 
conductivity. Grens [16] addresses the importance of these and other assumptions for the 
prediction of steady-state current density with a one-dimensional model. For the system 
to be considered here (i.e. a simple redox reaction in a well-supported electrolyte) it is 
reasonable to neglect these effects. 

The assumption of equal diffusivities (assumption 6) eliminates the diffusion equa- 
tion for one of the redox species. In a previous study [5] we included the effect of differing 
diffusion coefficients in the calculations for reversible LSV on a planar electrode, and 
found for a cathodic sweep when the diffusion coefficient of the oxidized species was 50% 
higher than that of the reduced species (a significant variation) the peak current was 
less than 2% greater and the peak potential 5/n mV (T = 298 K) more negative than 
if the diffusion coefficients were equal. If the kinetics are irreversible, the concentration 
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distribution of only the reactant species is important and assumption 6 has no effect 
on the voltammogram. Neglecting uncompensated resistance (assumption 7) limits the 
usefulness of the results if the potential drop between the reference electrode and the 
pore aperture is appreciable since the sweep rate at the pore opening would no longer 
be linear. However, uncompensated resistance may be negligible even when ohmic losses 
inside the pore are appreciable if the cross-sectional area available for current flow is 
larger outside the pore then inside, and the distance from the reference electrode to the 
pore mouth is small. Neglecting double-layer charging (assumption 8) places a lower 
limit on the reactant concentration since only the reaction current decreases with a de- 
creasing concentration, and an upper limit on the sweep rate since the charging current 
is a stronger function of the sweep rate than is the reaction current. 

Assumption 9 is an attempt to account for mass-transfer resistance in the bulk 
electrolyte by conceptualizing the pore aperture as an electroactive disk embedded on 
an insulating plane. Bond et al. [17] calculated the Bteady-state current to such a 
disk when the reactant surface concentration was uniform, and using Farady’s law they 
related the concentration gradient at the surface to the difference between the bulk and 
surface concentrations. Assumption 10 is invoked for mathematically convenience; in 
reality the potential at the tube inlet is radially dependent. A rigorous treatment which 
would eliminate the last two assumptions would require the simultaneous solution of the 
governing equations inside and outside the tube — a fairly intractable problem. 
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2.2 Governing Equations 


The concentration of the oxidized species is governed by the time-dependent, two- 
dimensional diffusion equation 


QCo _ PCo 1 d ( dCoX 

dt [ dz* + r dr \ dr ) 


( 2 ) 


where Co is the oxidized species concentration, V is the diffusion coefficient, t is the 
time, and r and z are the radial and axial coordinates, respectively. 

The concentration of the reacting species is initially uniform throughout the pore 
and equal to that in the reservoir (equation 2a). The material flux at the pore inlet is 
proportional to the difference between the bulk and inlet concentration (equation 2b) 
with the proportionality constant obtained by modeling the pore aperture as an electro- 
active disk on an insulating plane [17]. The concentration gradients at the back of the 
pore (equation 2c) and along the centerline (equation 2d) are zero, and the gradient at 
the pore wall is proportional to the reaction current (equation 2e). These constraints 
can be expressed mathematically as 


2a) 

t = 0; 

II 

£ 



2b) 

z — 0; 

= (-s/«0 (08 - Co) 

2d) 

o 

II 

CO 

o 

II 

*- 

2c) 

z = L\ 


2e) 

ii 

£ 

cT 

II 


where the superscript o indicates the spatially uniform initial concentration, F is Fara- 
day’s constant, and t n is the local reaction current density at the electrode surface which 
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is positive for a cathodic reaction. For reactions in which the kinetics are governed by 
the Butler- Volmer rate equation, *„ is given by 

*„ = nFfc 0 {Co exp [-anf(E w - E*)] -Cijexp [(1 - a)nf(E v - £“')] } (3) 

where k° is the standard rate constant, / = a is the cathodic transfer coefficient, 
Ey, and E°' are the solution potential at the tube wall and the formal potential of the 
redox couple, respectively, and Cr is the reduced species concentration. 

Due to the assumption of equal diffusivities, the concentration of the oxidized and 
reduced species are related by 

c o + Cr ~ Co(z,r,t) + C R (z,r,t) (4) 


The Laplace equation governs the solution potential 


d 2 E 

0 = e? + 


ld_( dE\ 
t 9r \ dr ) 


( 5 ) 


The potential at the pore opening starts at the equilibrium potential, E and changes 
linearly with time as set by the sweep rate, v = — ^ (equation 5a), the potential 
gradients at the closed end of the pore (equation 5b) and the centerline (equation 5c) 
are zero, and the gradient at the electrode surface is proportional to the reaction current 
(equation 5d). These boundary conditions can be expressed mathematically as 


5a) z = 0; (E — E^) = -vt 
5b) z = L\ f = 0 


5 e) r = 0; ^ = 0 

_ . . . QE , 

5d) r - d/2 ; k— = t n 
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where k is the conductivity of the electrolyte. 


Equation 2 is rewritten in a dimensionless form as 

dC* 0 PCI 1 6 ( Y dC b\ 

dr dX> Y dY V &Y ) 


6a) 

r = 0; 

C* 0 = 1 



6b) 

X = 0; 

^ = ( — 4/7r) (1 - C'o ) 

6d) Y = 0; 

**=0 

dY 

6c) 

x = r , 

**.0 

dX 

6e) y = l; 

dC '° - rzr 

dY v* 1 


where the dimensionless parameters are defined in Table 1, and I* is the dimensionless 
local reaction current density at the electrode surface 


7* = 


nFC° 0 {nfvVyi* 


( 7 ) 


The dimensionless form of equation 5 is 


d 2 E* 1 d ( y dE*\ 
~ dX* + YdY\ Y dYJ 
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8d) 


dE* 
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dE* 

dY 


= 0 

= (J/A)7* 


( 8 ) 


The dimensionless Butler-Volmer equation is 


I* = A 


C*ne~ aE 




( 9 ) 
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7 

ohmic resistance 

n 2 fFk°dCl 

J 

kinetic resistance 

2k 

a 

cathodic transfer coefficient 


7 

length 

radius 

21 

d 

A 

mass-transfer resistance 

k° 


kinetic resistance ( nfuV ) 1 ^ 2 

£ initial concentration ratio Cq/C% 


a 


diffusion time (Pnfv 

sweep time 4P 


Table 1: The definitions of the dimensionless parameters that govern LSV in a cylindrical 
electrode and their physical interpretation. 
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Since LSV in a cylindrical electrode depends upon all six parameters listed in Ta- 
ble 1, and each finite-element numerical simulation takes on the order of an hour of CPU 
time on an IBM/3090, it is desirable to make some simplifications to the above model, 
as discussed below. 

2.3 Negligible Axial Diffusion: With and Without Appreciable 
Solution Resistance 


The first simplification is to neglect axial diffusion. This has practical significance if the 
duration of a single linear Bweep is short since radial diffusion will predominate for the 
major fraction of the pore volume. Equation 6 reduces to 

dC'r 


dCb = 1 d (y dC'p\ 
dr Y dY\ 8Y ) 


( 10 ) 


10a) r = 0; C' 0 = 1 


10b) Y = 0; 


dCb 

BY 


= 0 


10c) y = i; -fir 


Taking the Laplace transform of equation 10 with respect to r, solving the resulting 
ordinary differential equation for the transformed concentration, relating the gradient 
of the transformed concentration at the electrode surface to the transformed reaction 
current, and finally applying the convolution property of Laplace transforms to invert 
back to the time domain results in 



( 11 ) 
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where A n , are the zeros of the first-order Bessel function Ji(A ni ) = 0. Equation 11 is 
valid regardless of the reaction rate equation since it relates the concentration of O at 
the electrode surface to an, as yet, unspecified reaction current. 


2.3.1 Reversible Kinetics 

For a reversible electrochemical reaction, the surface concentrations of O and R tire 
coupled through the Nernst equation to the solution potential by 

= (12) 
where E°' is the formal potential. Writing equation 12 in dimensionless form and utilizing 


equation 4 gives 


C‘o( l.r) = 


(1 + 


(13) 


{[! + .«] 

where E £ is an arbitrary function of time. Substituting equation 11 into equation 13, 
and solving for the integral expression results in a Volterra integral equation of the first 


kind 


f [l + Ee-^-^ 

Jo L n=! 


0 T ') <iT ' = 


t _ e El 

* srr =9(K) 


2(^[1 + e*i] 


(14) 


where the subscript rev denotes reversible kinetics. 


The dimensionless current, J* ct (t), can be determined numerically [2] by approxi- 


mating the integral in equation 14 as 




Q-l 


1^(1-')^' « X) B M<i h ~ 

**=o 


(15) 
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where h is the width of each interval such that r = qh and r' = fih. The coefficients 
are determined by approximating the unknown function by linear expressions 

in the individual intervals, integrating the left side of equation 15, and solving for the 
coefficients. For fi = 0 this procedure results in 




and for fi from 1 to q — 1 


B„ = h + £ 


..i *1 fc 

n=l nj 


(i _ „-«,<■) 1 


Substituting equation 15 into equation 14, and solving for I(qh) gives 


A-i 


g(K) - E - M) 


IUlh) = 


M=1 


Bo 


(16) 


(17) 


(18) 


If it is further assumed that the solution resistance is negligible, then the potential 
at the tube wall is equal to the potential applied at the pore opening (J5£ = ln£ - <tt), 
and the peak current and potential is obtained by stepping equation 18 through time 
until a maximum in the current is obtained. The voltammograms generated by the 
numerical solution of equation 18 depend only upon a and 

In order to assess the accuracy of the integral approximation given in equation 15, 
an alternative solution to equation 10 was generated by solving the equation for a unit 
step change in the potential at the tube wall and using Duhamel’s superposition principle 
[18] to account for the time- dependent boundary condition given by equation 13. The 
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resulting current-potential expression is 


C.M = 2(1 + f ) V? E _( 

n=l 


r exp 


- W T ' + *n.( T ~ T ')}] ^ 


1 + ^exp(-<rr , )] J 


(19) 


where A no are the zeros of the zeroth-order Bessel function J 0 (A n# ) = 0. 


In the presence of appreciable solution resistance, the two-dimensional Laplace equa- 
tion for the potential (equation 8) must be solved to obtain E £. Winsel [11] and Viner 
and Fedkiw [12] solved equation 8 to determine as a function of axial position for an 
arbitrary time- dependent potential variation at the pore opening. For a linearly varying 
potential at the pore mouth, the result is 


El(X,r) = sin/Ur f /*(*', r) sin /UT'dJT (20) 

A 7 „=l Pnh\Pn) JO 

where 1 0 (/3 n ) and 1 1 ( /3 n ) are zero- and first-order modified Bessel functions, respectively, 
and j3 n = (n — \)n fa- Equation 20 is valid regardless of the reaction rate equation. 


2.3.2 Irreversible KineticB 

When the kinetics are irreversible, the back reaction may be ignored and for a cathodic 
reaction the Butler- Volmer kinetic expression given in equation 9 reduces to 

= A (21) 


98 




Correlations which relate the peak currents and potentials obtained from LS\ in a 
tubular electrode to the governing dimensionless parameters are obtained for reversible 
(Nernstian) and irreversible kinetics by assuming solution resistance and axial diffusion 
are negligible (equations 18 or 22). Solution resistance is then accounted for by coupling 
the analytical solution to the one-dimensional diffusion equation to the analytical ex- 
pression for the solution potential at the electrode surface (equation 20). In parameter 
regions in which numerical difficulties prevent these coupled non-linear equations from 
being evaluated, a finite-element numerical solution is used to solve equations 6 and 8. 
The parameter regions in which axial diffusion may be ignored is determined by com- 
paring the peak currents obtained with and without the axial component of the diffusion 
equation included in the model. 
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3.1 Negligible Axial Diffusion 


3.1.1 Negligible Solution Resistance 
Reversible kinetics 

In the absence of solution resistance, equation 18 was applied by letting Eh = 
(ln£ — err) for 51 different values of <r ranging from 10 -1 to 10 4 at 30 different values 
of £ ranging from 10 -3 to 10 s . Richardson extrapolation [20] was used to speed the 
convergence of the infinite series for the coefficients (equations 16 and 17), and all 
other non-oscillating series given in this paper. Figure 1 shows the dimensionless peak 
currents and potentials as a function of y/a for three different values of £. Only some of 
the calculations up to y/cr = 20 are shown in Figure 1 in order not to clutter the figures 
and to emphasize the behavior at large and small a. As a — > oo, the diffusion layer is 
small relative to the pore diameter and the peak values approach those obtained from a 
planar electrode (i.e. the peak current is proportional to v 1 !* and the peak potential is 
sweep rate independent). Weidner and Fedkiw [5] developed empirical correlations for 
the peak current and potential as a function of £ for LSV on a planar electrode and these 
are listed in Table 2. The peak current obtained from LSV on a planar electrode is 
within 5% of the peak current obtained from that on a cylindrical electrode for a > 300, 
and the peak potentials are within 5/n mV of each other for a > 40. 

As the sweep rate approaches zero (<t —* 0), the radial concentration gradients 
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0 5 10 15 20 


Figure 1: Dimensionless peak currents and potentials for linear-sweep volta mm etry in 
a cylindrical electrode when the kinetics are reversible, and ohmic losses and the axial 
diffusion flux are negligible. The symbols result from the use of equation 18 and the solid 
lines are the correlations which were fit to the calculations (equations 24 and 26). The 
dashed lines as<r-»oo are the dimensionless peak currents and potentials obtained from 
LS V on a planar electrode [5] , and those as u — > 0 are the equations for the dimensionless 
peak currents and potentials in a thin-layer cell [6] (see Table 2). The potential, in mV, 
i6 for T = 298 K. 
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(A m ) 



Peak Current or Potential Expression 


Reference 


<r Limit 


£ > i; 

^Vev(thin-layer) = 

*<i; 


8* 

2(1 + 0 


[ 6 ] 


a < exp 3.91 


( oWTt )- 2 - 81 


e > i ; 0 

^,r e v( ttin - la yer) = 

£<i; i»£ 


[ 6 ] 


<r < 0.5 


■^p r e« (planar) = 0.6103-0.164 


(l.os+e) 


[5] 


<r > 300 


(P lanar ) = In 




1.19 


I 0.84 


+ 3.74 


[5] 


<r > 40 


Table 2: Expressions for the peak current and potential on a planar electrode and in a 
thin-layer cell for reversible kinetics. The er limit indicates the value at which the peak 
currents are within 5% and the peak potentials within 5/n mV (T = 298 K) of the values 
calculated from LSV in a cylindrical pore if solution resistance and axial diffusion are 
negligible. 
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approach zero and the voltammograms have the characteristics of LSV in a thin-layer 
cell. Hubbard and Anson [6] derived a current-potential relationship for LSV in a thin- 
layer cell when only the oxidized species is initially present (£ = oo) and the voltage is 
swept in the negative direction. They found that the cathodic current is proportional to 
the sweep rate and given by 


I rev (thin-layer) = 


n 3 f¥vVC$e E ' 

A [1 + «*•]' 


(23) 


where V is the volume of the thin-layer cell. When the initial concentration of the reduced 
species is not negligibly smalj, equation 23 remains valid except that Cq is replaced by 
(i Cq + Cr). When Cq > Cr (£ > 1), the peak is obtained by differentiating equation 23 
with respect to time and setting the resulting equation equal to zero. The potential 
which satisfies the resulting equation is the formal potential ( E * = 0) which, when 
substituted into equation 23, gives the peak current. However, when Cq <Cr((< 1) 
the starting potential is negative with respect to the formal potential and therefore the 
peak current occurs at the start of the sweep ( E * = ln£) and the current decays as the 
oxidized species is reduced. We see from Figure 1 that the dimensional peak current 
is proportional to the sweep rate at low <r and the peak potential approaches zero as 
predicted from equation 23. The dimensionless peak values in the limit of <r — * 0 are 
listed in Table 2. The upper limit of a for which the the thin-layer-cell approximation 
is valid is a strong function of £. For ( > 5, the peak current obtained from equation 18 
is within 5% of Ip <rev (thin-layer) given in Table 2 when <r < 2, but for £ < 0.1, <j has to 


103 



» 

be less than 0.1 for the two currents to be within 5% of each other. A one-parameter 
least-squares fit was performed for the upper <t limit as a function of *, and the result is 
shown in Table 2 


The following empirical relationship, motivated by the asymptotic trends seen in 
Figure 1, was fit to the dimensionless peak currents calculated from equation 18 

= [^,rev(tHin-layer) ro + j; iW (planar) m ] 1/m (24) 

where m is an empirical parameter which was found by least-squares fit of equation 24 
to the dimensionless peak currents for 51 different values of <r. The parameter m was 
obtained at 30 different values of *, and a two-parameter empirical equation for m was 
generated as a function of * and is listed in row 1 of Table 3. 


Substituting I* rev (thin-layer) for * > 1 and /* >rTO (planar) into equation 24 along 
with the definitions for I* and cr gives 


(*)'- 


wn*fF<PLC&{ 1+01 
16 * 


v™!' + [l.57(nF) 2 d 3 i(7S/; iret ,(planar)] m (25) 


where t,, is the peak current which would be measured if solution resistance was negligible. 
Therefore, plotting {i v / versus j/ TO / a should yield a straight line if the kinetics are 
reversible, and the slope and intercept may be used to estimate n. Table 3 lists the 
proper way to plot LSV data in order to yield a straight line, along with the equations 
for the resulting slope and intercept. 


The following empirical relationship, also motivated by the asymptotic trends seen 
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in Figure 1, was fit to the theoretically predicted peak potentials from equation 18 


E; tTt v(<r) = [E;, fev (: planar) - J?; >rev ( thin-layer)] [^-^] + E; >r „ (thin-layer) (26) 

The two- parameter empirical equation for m as a function of £ is listed in row 2 of 
Table 3. 

Only one parameter (m) was used to fit the peak current and potential calculations 
in order to facilitate a graphical analysis of experimental data. The disadvantage to 
a one-parameter fit is that equation 24 is only accurate to within 5% and equation 26 
to within 5/n mV (T = 298K). Although this error may be significant, peak values 
estimated from equations for LSV on a planar electrode or in a thin-layer cell result in 
even larger errors if the value of a is outside the range listed in Table 2. Therefore, these 
correlations provide more accurate estimates of the peak values at moderate sweep rates 
than is possible with planar or thin-layer theory alone. 

The error in the simulated voltammogram arising from the integral approxima- 
tion given in equation 15 was assessed by comparing the peak currents obtained from 
equation 18 to those obtained from equation 19 which was derived using Duhamel’s su- 
perposition principle. Although the peak current for equation 19 may be obtained by 
generating an entire voltammogram, a more convenient analytical procedure was derived 
by differentiating equation 19 with respect to r and setting the result equal to zero. The 
potential which satisfies the resulting expression is the peak potential which, when sub- 
stituted into equation 19, specifies the peak current. Newton’s method was used to solve 
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the resulting non-linear equation for the peak potential generated from this procedure. 
At ( = 1000, the peak currents for cr varying from 10 -1 to 10 4 calculated from the two 
current-potential relationships agreed to within 0.1% of each other, verifying the integral 
approximation . 

Irreversible kinetics 

In the absence of solution resistance, equation 18 was applied by letting E £ = 
(ln£ — <tt) for 51 different values of u ranging from 10 -1 to 10 4 at 9 different values of 
a ranging from 0.1 to 0.9. As stated by Matsuda and Ayabe [3] and reaffirmed in our 
calculations, the kinetics can be assumed irreversible if A < 10 -s and was the value used 
in the calculations. Figure 2 shows the dimensionless peak currents and potentials as 
a function of \fa for three different values of a. As with reversible kinetics the peak 
values approach those obtained from a planar electrode as a — ► oo [4] and approach 
the thin-layer peak values as <r — ► 0 [7]. Table 4 lists the peak current and potential 
expressions for LSV on a planar electrode and in a thin-layer cell which are used in the 
correlations. 

The following empirical relationship, which is analogous to that used when the 
kinetics were reversible, was fit to the theoretically predicted dimensionless peak currents 
obtained from equation 22 

Jp.irrevH = fcrrev (ton-layer)"* + i; (irr „(planar) m ] ^ (27) 
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Figure 2: Dimensionless peak currents and potentials for linear-sweep voltammetry in 
a cylindrical electrode when the kinetics are irreversible and ohmic resistance and axial 
diffusion are negligible. The symbols result from the use of equation 22 and the solid lines 
are the correlations which were fit to the calculations (equations 27 and 28). The dashed 
lines as a -> oo are the dimensionless peak currents and potentials obtained from LSV 
on a planar electrode [4], and those as a -> 0 are the equations for the dimensionless 
peak current in a thin-layer cell [7] (see Table 4). The potential ordinate approaches 
(— In <r) as <r — ► 0. The potential, in mV, is for T — 298 K. 
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Table 4: Expressions for the peak current and potential on a planar electrode and in 
a thin-layer cell for irreversible kinetics. The <r limit indicates the value at which the 
peak currents are within 5% and the peak potentials within 5/n mV (T = 298 K) of the 
values calculated from LSV in a cylindrical pore if solution resistance and axial diffusion 
are negligible. 
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were m is an empirical parameter which was found by least-squares fit of equation 27 to 
the dimensionless peak currents for 51 different values of <r. The fit was performed at 
nine different values of a, and all the m values were within 1% of each other. Taking m 
as a constant, a can be estimated from the slope of a versus plot. The 

Blope and intercept of the straight line and the average value for the parameter m are 
listed in row 3 of Table 3. 

The following empirical relationship, motivated by the asymptotic trends seen in 
Figure 2, was fit to the theoretically predicted dimensionless peak potentials obtained 
from equation 22 

exp aE*i rrev ((r) = [(exp aE; >irrcv (thin-layer)) " + (expa£?; >tVr „(planar)) m ] 1 (28) 

where m is an empirical parameter which was found by least-squares fit of the above 
equation to the peak potentials. The nine m values, corresponding to the different a 
values, were within 5% of each other with the average value listed in row 4 of Table 3. 
Equation 28 may be simplified and the peak potentials written in the form which is 
plotted in Figure 2 giving 

aE; iirrev - InfA/v'S) = 0.488 In 0.202 + ==) j (29) 

After a is determined from the experimental peak-current data, the peak potentials can 
be used to estimate the standard rate constant, k°, by plotting the data as indicated in 
row 4 of Table 3 and evaluating the slope and intercept. 
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3.1.2 Appreciable Solution Resistance 
Reversible kinetics 

Since the current and the potential are position dependent when solution resistance 
is appreciable, and the solution potential at the tube wall (equation 20) is a function of 
the total current, equations 18 and 20 are coupled and must be solved simultaneously. 
The total current as a function of applied potential at each time step was obtained 
by dividing the pore into N segments of equal length and calculating the current in 
each segment using equation 18 with the total current being the sum of the individual 
segment currents. The N current and N potential equations were simultaneously solved 
using the IMSL [19] non-linear equation-solver subroutine DNEQNJ which makes use of 
a user-supplied Jacobian. The Filon algorithm [21] was used to perform the integration 
in equation 20, and a binomial averaging algorithm [22] was used to extrapolate the 
infinite series. Unfortunately, numerical difficulties prevented the non-linear equation 
solver from converging in a parameter region were solution resistance was appreciable. 
Increasing the number of divisions in the pore may have allowed the non-linear equation 
solver to converge, but at N = 20 the amount of computer time per simulation was 
on the order of that needed to solve numerically the two coupled partial differential 
equations, and so the following procedure was used to account for solution resistance 
when the kinetics are reversible. 
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Equations 6 and 8 were simultaneously solved using the IMSL [19] finite-element 
procedure PDE/PROTRAN, and the Butler- Volmer kinetic expression (equation 9) was 
numerically integrated along the tube wall using the IMSL [19] subroutine DQDAGS at 
each time step. The Butler- Volmer equation is valid for reversible kinetics as long as 
A is large. Matsuda and Ayabe found that the kineticB could be assumed reversible if 
A > 15, and so a value of A = 50 was used in the calculations to insure reversibility. It 
was observed that for a = 1 and £ = 100 the peak currents and potentials depend only 
upon the grouping (J/ A)j 3 = 0 as long as 7 was greater than 100. The lower bound 
on 7 is due to appreciable axial diffusion of material into the pore through the pore 
aperture. The effect of axial diffusion on the voltammogram will be addressed further 
in the next section. 

The peak currents, normalized by the peak currents which would be obtained if 
solution resistance was negligible, axe plotted in Figure 3 as a function of 0. From 
Figure 3 we see that at 0 = 1 solution resistance diminishes the peak current by less 
than 1%, but at 0 = 10 the peak current has been reduce by 13% due to ohmic effects, 
and by 24% at 0 = 20. 

In order to assess the numerical accuracy of the finite-element numerical solution, 
an alternative current calculation was performed by numerically integrating the current 
in the electrolyte, i t , at the inlet to the tube using Ohm’s law (i t = kVE). The potential 
gradient was estimated by using a three-point-derivative formula. The two independent 
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Figure 3: The normalized peak current in a cylindrical electrode when the kinetics 
are reversible as a function of the dimensionless grouping (J/A)7 = The peak 

current, t p , is normalized with respect to the peak current which would be obtained if 
solution resistance was negligible, i p . The symbols are the results from the finite-element 
numerical solution, and the Bolid line is the correlation which was established for a 
thin-layer cell [24] when <r = 1. 
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current calculations agreed to within 2% of each other for the results given in Figure 3, 
but for 0 > 20 the the discrepancy between the currents was larger than 2%. 

As shown in Section 3.1.1, the voltammogram will have characteristics of LSV in a 
thin-layer cell when cr < 3 in the absence of solution resistance. Hinman et al. [23] used 
a network model to simulate LSV in a thin-layer cell to account for ohmic distortion 
when the kinetics are reversible. Three assumptions, in addition to the ones made in 
obtaining the results shown in Figure 3, were made in their analysis: (1) the potential 
is radially uniform, (2) only the oxidized species is present initially (£ = oo), and (3) 
the rate at which the potential changes with time is uniform throughout the cell as set 
by the sweep rate. Using the same set of assumptions we applied a continuum model to 
the thin-layer cell [24], and found that the normalized peak current was a function of a 
single variable which in the present nomenclature is 0^/ir. The resulting correlation is 
plotted in Figure 3 for <r = 1. 

Both models predict a precipitous drop in the peak current as 0 is increased past 
one, but the drop is steeper in the present calculations. For example, at 0 = 10 the peak 
current calculated from the thin-layer model [24] is 9% higher than that calculated by 
numerically solving the two coupled partial differential equations (equations 6 and 8). 
The assumption of a radially-uniform potential should be valid at <r = 1 since the thin- 
layer approximation holds for a < 3 (see Tables 2) and a small concentration gradient 
leads to a small potential gradient. The difference in the normalized peak currents 
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resulting from a non-zero initial concentration of the reduced species should be negligible 
since the difference between the peak current for £ = 100 and £ = oo in a thin-layer cell 
when solution resistance is negligible is only about 1%. The difference in the peak 
currents predicted from the two models is therefore expected to result from assuming 
the sweep rate is uniform throughout the pore. As indicated by equation 23, the local 
current is not only a function of the axial dependent potential, but it is also directly 
proportional to the local rate of change of the potential, d\E\/dt. Since d\E\/dt decreases 
with distance into the pore, setting d\E\/dt equal to the sweep rate, u, overestimates the 
current. 

Yet to be determined is the functional dependence of the normalized peak current 
on o and 0. From our previous study [24] we expect the normalized peak currents to 
only be a function of the grouping Qy/o when the thin-layer approximation holds. The 
shift in the peak potential due to ohmic effects also needs to be quantified. The effect of 
£ on the peak current in the presence of solution resistance has not been explored either, 
but the effect on the normalized peak current is expected to be small. 

Irreversible kinetics 

The effect of solution resistance on the voltammogram when the kinetics are irre- 
versible and axial diffusion is negligible was determined by simultaneously solving equa- 
tions 20 and 22 with the pore divided into 10 segments. This value for N was used since 
no noticeable effect on the peak current was observed for larger N, and the computer 
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time per simulation scaled as AT 2 . The normalized peak currents were found to depend 
only upon 0 for <r > 20, and Qy/a for cr < 2, Figures 4 and 5 show the normalized peak 
currents obtained for <r > 20 and <r < 2, respectively. For a > 20, numerical difficulties 
prevented the non-linear equation solver from converging for 0 > 10, and so in this © 
range equations 6 and 8 had to be solved using the finite-element numerical solution. 
Comparisons between the two solution procedures were made at 0 < 10 for a = 100 
and the resulting peak currents were within 0.2% of each other as long as 7 > 10, a 
difference which would be visually indistinguishable in Figure 4. All the calculations 
shown in Figure 5 were obtained from the finite-element solution. 

A totally empirical, three parameter least-squares-fit of the results in Figure 4 pro- 
duced the following correlation for the normalized peak current as a function 0 for 
a > 20 . 

(tp/tp).Vrev = [(1.372 - 0.438 log 0)- 7 - 50 + 1] -0 ' 1333 (30) 

Five different values of a between 20 and 500 were used in the correlation in which 
0 r an ged from 10 -2 to 10 2 . A three parameter least-squares-fit was also made of the 
calculations in Figure 5 which resulted in the the following correlation for the normalized 
peak current as a function ©■y/tr for a < 2. 

(»p/»p)i™ = 1(1.606 -0.4971og©v^) ' +lj (31) 

In establishing this correlation, a ranged from 0.2 to 2, and 0 from 10 to 100. 

The simulations which established equations 30 and 31 were run with a = 0.5. It is 
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Figure 4: The normalized peak current in a cylindrical electrode as a function 0 for 
<r > 20 when the kinetics are irreversible and axial diffusion is negligible. The peak 
current, » p , is normalized with respect to the peak current which would be obtained if 
solution resistance was negligible, t p . The circles (o) are the results from simultaneously 
solving equations 20 and 22, and the diamonds (O) are the results from the finite-element 
solution. The solid line is the correlation which was fit to the calculations (equation 30). 
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Figure 5: The normalized peak current in a cylindrical electrode as a function 
for a < 2 when the kinetics are irreversible and axial diffusion is negligible. The peak 
current, t p , is normalized with respect to the peak current which would be obtained if so- 
lution resistance was negligible, i p . The symbols are the results from the finite-element 
solution, and the solid line is the correlation which was fit to the calculations (equa- 
tion 31). 
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not expect that normalized peak currents will be a strong function of a (if at all), but 
that remains to be determined. Analogous correlations to those for the normalized peak 
currents as a function of 0 still need to be established for the shift in the peak potential 
due to solution resistance. Also, more simulations are required in order to quantify the 
ohmic effects in the cr region between 2 and 20. 

In order to check the accuracy of the finite-element calculations, the currents ob- 
tained from the integration of the Butler- Volmer kinetic expression were compared to 
those resulting from radially-integrating Ohm’s law at the pore aperture, and for a > 20 
and 0 < 100 an agreement within 1% was obtained. No calculations were performed 
at larger values of 0 since at 0 = 100 the voltammogram is severly distorted with the 
peak spread over tens of mV, and bo simulating LSV beyond this point would be of no 
practical use. 

3.2 Effect of Axial Diffusion 

The above correlations will not hold if appreciable material diffuses into the pore from 
the reservoir which will occur if u is too small for a given 7. In order to determine the 
lower bound of <r as a function of 7 above which axial diffusion can be neglected, the 
finite-element technique was used to solve equation 6 with the potential at the pore wall 
set by the sweep rate. For the case of reversible kinetics, the Nernst equation was UBed as 
the boundary condition at the electrode surface instead of the Bulter- Volmer equation, 
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and it was observed that for £ = 10, the current due to diffusion of material into the 
pore contributes less than 2% to the peak current if cr > 42/7. This bound of a is valid 
for £ > 10, and results in a percentage error less than two for ( < 10. When the kinetics 
were irreversible, it was observed that axial diffusion can be ignored as long as c > 75/7- 

Material diffusing through the pore mouth will contribute a larger fraction of the 
current as the potential becomes more non-uniform since the potential driving force at 
the pore aperture will need to be increased in order to drive the reaction at the back of 
the pore to an appreciable rate. Therefore, a larger value of <r will be required at a given 
7 so that the assumption of negligible axial diffusion is valid. The relationship between 
<r and 7 has yet to be determined when solution resistance is appreciable. 

4 Summary 

Correlations were established for determining kinetic parameters from LSV in a cylin- 
drical electrode. Due to the complex dependence of the peak currents and potentials on 
the governing parameters, simplifications were introduced at various parameter limits. 
The first simplification was to neglect axial diffusion which allows the diffusion equa- 
tion to be analytically solved giving the surface concentration of the redox species in 
terms of the reaction current and the solution potential at the electrode surface. Fur- 
ther assuming solution resistance is negligible reduces the computation time over three 
orders-of-magnitude from the finite-element solution numerical procedure facilitated the 


120 



establishment of the desired correlations. In the limit of reversible (Nemstian) and 
irreversible kinetics, correlations were established which relate the peak currents and 
potentials to the sweep rate by making use of the peak expressions on a planar electrode 
and in a thin-layer cell. 

The correlations, which are established assuming solution resistance is negligible, 
are also useful when solution resistance is appreciable since the currents can be adjusted 
to give peak data which would have been obtained in the absence of ohmic effects. The 
effect of solution resistance on the peak current was quantified for the irreversible-kinetic 
case for <r > 20 and a < 2, but more simulations are needed to quantify the peak currents 
for c falling between these two limits. More calculations also need to be performed in 
order to establish correlations when the kinetics are reversible. Correlations have yet to 
be established for the peak potential calculations which has been collected thus far. 

Acknowledgement — This work was supported by a grant from the NASA-Lewis Research 
Center (NAG 3-649). 
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Notation 


electrode area, cm 2 
coefficient in integral approximation 

concentration of oxidized species of redox couple, mol/cm s 
concentration of reduced species of redox couple, mol/cm 3 
initial concentration of oxidized species of redox couple, mol/cm 3 
initial concentration of reduced species of redox couple, mol/cm 3 
Co/Cp 

diffusion coefficient of the redox species, cm 2 /s 
solution potential, V 
formal potential of redox couple, V 
equilibrium initial potential, V 
dimensionless potential, nf(E — E°') 

F/RT V/equiv 
Faraday’s constant, C/equiv 
width of interval in integral approximation 
total reaction current, A 

total reaction current when solution resistance in negligible, A 
current density at electrode Burface, A/cm 2 
current density in solution, A/cm 2 
exchange current density, A/cm 2 
zeroth-order modified Bessel function 
first-order modified Bessel function 
dimensionless current, * n / (nfCp y/ nfvV) 
zeroth-order Bessel function 
standard rate constant, cm/s 
empirical parameter 
number of e~ in redox couple 
gas constant, J/(mol-K) 
time, s 

temperature, K 

distance from center of pore, cm 
volume of pore, cm 3 
dimensionless axial coordinate, 2 zjd 
dimensionless radial coordinate, 2 r/d 
z distance from the pore aperture, cm 
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Greek 


A. 

^ni 

K 

V 

T 


(» - \)*h 

zeros of zeroth-order Bessel function J 0 
zeros of first-order Bessel function Ji 
conductivity of the electrolyte, (ohm • cm) -1 
sweep rate, — V/s 
dimensionless time, td/2V 


Subscripts 

irrev irreversible kinetics 

p value at peak of volt ammogr am 

rev reversible kinetics 

w value at tube wall 
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SECTION 6 


ON-GOING WORK 

The PhD candidate, Mr. Weidner, is currently performing linear sweep voltammetry 
experiments for a number of redox couples in a tubular electrode to confirm the theoretical 
calculations presented in Section 5. Upon completion of this work he will then repeat these 
same measurements for a nickel electrode and attempt to verify the methodology for this 
electrode reaction. We anticipate completion of the dissertation no later than May, 1991 
and we shall forward, when available, the dissertation and any additional papers which 
result from the project. 
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